1 Introduction
1.1 An Overview
A maritime transportation system is greatly impacted by the efficiency of the routing and scheduling approach being adopted. The United States, which is the largest trading nation in the world in terms of imports and exports, accounts for nearly 20% of the global trade, of which over 99% is comprised of sea cargo (Agarwal, 2007). U.S. ports and waterways handle about 2.5 billion tons of trade annually, and this is expected to double within the next 15 years (Agarwal, 2007; Al Khayyal and Hwang, 2007). The substantial worldwide increase in sea-borne trade has led to a proportional growth of the international fleet. Maritime transportation constitutes the majority of long-distance shipments in terms of volume (Rodrigue et al., 2017). Oil is the world’s leading transported product in terms of tonnage, amounting to about 4.8% of the global transported value. In 2015, the total exported value of crude oil was estimated at $786.3 billion, and the share of the Gulf countries was valued at $325 billion or 41.3% of the global crude oil export (Trade Map, 2017). The operation of a typical fleet of oil tankers is highly expensive due to the high costs of owning and chartering oil tankers in addition to the costs associated with the daily vessels’ operation. In recent years, the prices of oil have dropped from $120 per barrel (2012) to a current price of about $69 per barrel (June, 2018), drastically affecting economies of major oil producing countries such as Saudi Arabia and Kuwait. As a result, efficient fleet utilization has become a critical issue for shipping companies, which, in turn, has prompted researchers to develop effective optimization-based decision support systems for the management and routing of ships (see, e.g. Christiansen et al., 2004, 2013; Christiansen and Fagerholt, 2009; Hoff et al., 2010; Pantuso et al., 2014). Ronen (1983, 1993) indicates that the maritime transportation mode has attracted the least research among all transportation modes, principally due to the complexity of the problem and the high degree of uncertainty involved in shipping operations. Nonetheless, with advances in technology, there has been a growing interest in this area of research as highlighted by Christiansen et al. (2004).
1.2 Problem Statement
This paper investigates a single source-single destination crude oil vessel scheduling transportation-inventory problem faced by the Kuwait Petroleum Corporation (KPC), among a host of other transportation scenarios. KPC aims to transport substantial quantities of crude oil from Shuaiba Port in Kuwait (source) to the storage facility at Dalian Xingang Port in China (destination). Although oil companies face many demand-transportation related crude oil scenarios, we focus in this paper on a single-source and a single-destination operational scenario that is typically encountered by many large-scale oil companies such as KPC and the Saudi Aramco. This is motivated by the fact that for a destination having considerable long-term demand requirements, KPC adopts a long-lasting practice of allocating a separate fleet of vessels to meet such demand requirements without interfacing with the demands at other destinations. For example, KPC aims to meet demand requirements agreed upon with the Chinese government in concert with major oil companies (such as China Petroleum and Chemical Corporation or China National Petroleum Corporation) to transport large quantities of crude oil from Shuaiba Port in Kuwait to the storage facility at Dalian Xingang Port using a dedicated fleet of vessels. (To ease readability, we will henceforth use the terms “product, source, and destination” in lieu of “crude oil, Shuaiba Port, and Dalian Xingang Port”, respectively.)
Typically KPC transports large quantities of crude oil, and so it is practically more advantageous for KPC to use fully loaded vessels at the source to be fully discharged at the destination. The daily export from the source depends on the availability of the product at the source, storage capacities at the destination, rates of consumption of the product at the destination, the availability of vessels, and the current storage level at the destination storage facility. KPC can also lease a transshipment facility, if such a decision enhances the efficiency and cost effectiveness of the overall operation. The demand at the destination is governed by the consumption rates at the destination, which might not be fixed during the entire time-horizon and might vary based on, for example, seasonal considerations. The level of the product at the destination storage facility is desired to lie within certain lower and upper bounds, and hence, a daily penalty is imposed based on limited permissible shortage or excess quantities with respect to these bounds. The fleet of vessels used by KPC for transporting the product from the source to the destination is composed of different types of vessels, some of which are self-owned while others can be possibly chartered for a specified duration of the time-horizon. Each vessel-type is characterized by its size, speed, loading and unloading times, etc. Accordingly, KPC aims to satisfy the demand requirements while minimizing the total cost that is comprised of the operational expenses of vessels, maintenance and leasing expenses of the transshipment facility, penalties resulting from violating desired lower and upper limits on the storage level, and the cost of chartering vessels.
1.3 Contribution and Organization
The present paper extends the work of the authors (Sherali and Al-Yakoob, 2006a, 2006b; Al-Yakoob and Sherali, 2013) related to vessel scheduling. Exact and aggregated mixed-integer programming formulations and related rolling horizon heuristics for the same demand structure considered herein were presented in Sherali and Al-Yakoob (2006a) for the single source-destination operation scenario, and in Sherali and Al-Yakoob (2006b) for the multiple sources and multiple destinations, while considering the leasing of temporary transshipment storage facilities, each of which has a known discretized location. It was shown that schedules obtained via the modelling approach adopted therein were superior to those generated by a manual approach in a case study related to KPC. Also, it was emphasized in Sherali and Al-Yakoob (2006b) that the use of transshipment facilities has the potential of reducing overall operational costs of the generated schedules, because this permits shipments to be stocked closer to the destinations during slack periods for use during periods of higher demand, thus resulting in using fewer vessels and, in particular, circumventing the high costs associated with chartering vessels. Subsequently, Al-Yakoob and Sherali (2013) solved the problem studied in the latter paper using a column generation approach, which led to further cost reductions. This paper examined the impact of allowing a transshipment facility to be situated at one of a number (three, specifically) of plausible discretized locations, instead of just a single specified location as done in Sherali and Al-Yakoob (2006b), and demonstrated that further cost reductions can be thus achieved. This observation motivated the idea that examining a continuum of potential locations for a given storage facility could possibly further reduce overall operational costs quite significantly. The present paper extends the foregoing body of work by including modelling considerations for determining a cost-effective location for the transshipment facility within a restricted continuous region. In this context, we design a hybrid two-stage algorithm that inputs a set of discrete judicious facility locations resulting from solving the proposed model using a specialized rolling horizon algorithm into the column generation method of Al-Yakoob and Sherali (2013).
The remainder of this paper is organized as follows. Section 2 reviews the related literature. We introduce some modelling preliminaries in Section 3. A mathematical programming model is presented in Section 4 for the single source-destination problem with a restricted continuous region for locating a single transshipment facility. Specific modelling considerations for the location aspect of the transshipment facility are discussed in Section 5. To further reduce overall operational costs, a hybrid, two-stage algorithmic approach is proposed in Section 6. Section 7 then provides computational results using a wide range of simulated test cases representing different operational scenarios related to the Kuwait Petroleum Corporation in order to assess the impact of the transshipment facility on the cost and efficiency of the overall operation. Finally, Section 8 concludes the paper with a summary and recommendations for implementation and future research.
2 Related Research
Several modelling and algorithmic approaches have been employed in the literature to solve vessel transportation problems, principally mathematical programming methods (see, for example, Persson and Göthe-Lundgren, 2005; Ronen, 2002; Sherali and Al-Yakoob, 2006a, 2006b; Al-Yakoob and Sherali, 2013). A particular issue that we aim to explore in the present paper relates to finding an optimal location for a transhipment facility within the framework of a vessel transportation-inventory operation, which results in a unique facility location-allocation plus routing and scheduling problem. The classical discrete or continuous location-allocation problem (Aykin and Brown, 1992; Sherali and Adams, 1984; Shetty and Sherali, 1977) seeks to determine the location of a number of supply centres in a plane to serve several customers having fixed locations, and to simultaneously specify the flow allocations of these supply centres, with the objective of minimizing the total distribution costs. Different variants of this problem have been extensively investigated in the literature and many modelling approaches and algorithms have been devised to solve such problems – we refer the reader to Drezner and Hamacher (2004) for a detailed survey of this fertile research area. However, none of the existing work addresses the location-allocation problem within the context of a vessel scheduling problem as addressed herein.
We next summarize some of the relatively more recent vessel scheduling papers published in the open literature. In Xinlian et al. (2000), the authors studied a fleet planning problem that seeks to specify types of ships to add to an existing fleet of ships along with determining an optimal fleet deployment plan. Fagerholt and Christiansen (2000a) investigated a multi-product scheduling problem, where each ship in the fleet is equipped with a flexible cargo hold that can be partitioned into multiple holds in a specified number of ways. In follow-on research, Fagerholt and Christiansen (2000b) employed a set partitioning approach for this problem and designed an algorithm for finding optimal schedules for individual ships. A mixed-integer programming model and a two-stage cost-based heuristic procedure was developed in Ronen (2002) for a vessel transportation problem. The problem of allocating bulk cargoes to tanks in maritime shipping was studied in Hvattum et al. (2009). A special mixed-integer programming model was formulated in Furman et al. (2011) for a transportation scenario faced by ExxonMobil to ship volumes of vacuum gas oil from origin points located in Europe to refineries located in the United States. Andersson et al. (2015) studied the impact of speed on the planning of shipping routes, and proposed an algorithm for a deployment and routing problem in roll-on/roll-off (Ro–Ro) shipping scenario, in which case the vessels can accommodate wheeled cargoes (cars, trucks and trains). Hennig et al. (2011) proposed a path flow framework with a priori route generation for a crude oil operational scenario accommodating splitting cargoes. Hennig et al. (2012) presented another path flow model for an oil tanker routing and scheduling problem with the objective of finding load sizes, ship routes as well as port arrival and departure times that minimize the overall transportation expenses. Hennig et al. (2015) investigated an oil tanker routing and scheduling scenario with split pickup and split delivery, and studied the applicability of the path flow-modelling methods (developed in Hennig et al., 2011, 2012) within a column generation framework.
Aizemberg et al. (2014) studied a petroleum transportation problem, and presented a column generation-based framework to solve the problem. Al-Yakoob and Sherali (2013) presented a column generation heuristic for a mixed-integer programming model to find an optimal combination of vessel schedules in a multiple sources and multiple destinations scheduling scenario. Kobayashi and Kubo (2010) presented a column generation algorithm for a time–space network to solve a vessel routing and scheduling problem, where cargoes can be assigned to multiple compartments onboard vessels. Li and Pang (2011) and Pang et al. (2011) investigated a shipping problem where each terminal/berth handles only one vessel at a time, and proposed a set partitioning formulations and column generation techniques to solve the problem. Halvorsen-Weare et al. (2012) presented a voyage-based approach for a supply vessel-planning problem.
A problem of transporting products from oil refineries to storage facilities was studied in Persson and Göthe-Lundgren (2005), where inventory levels of both refineries and facilities are affected by the process of scheduling at the refineries and by the demand at the facilities. A column-generation modelling framework was designed to solve the problem. Also, a Dantzig–Wolfe column generation approach was employed in Brønmo et al. (2010) to solve a multi-ship pickup and delivery problem with time-windows and flexible cargo sizes. A constraint programming-based heuristic procedure embedded within a column generation framework was proposed in Pang and Li (2011) for a ship routing problem while considering loading and unloading times of cargoes at pickup and drop-off locations. A flexible modelling approach for the inventory routing problem was presented in Song and Furman (2013), which can accommodate different practical issues. The authors proposed an optimization based heuristic method to solve the problem. Other related inventory routing problems have been investigated in Al Khayyal and Hwang (2007) and Christiansen and Fagerholt (2009). Soroush and Al-Yakoob (2018) extended the work of Sherali and Al-Yakoob (2006a) by incorporating stochastic aspects into the modelling approach presented therein for the same vessel transportation-inventory scenario studied in the current paper.
Finally, truck-petrol-station distribution planning and replenishment-related problems have been extensively investigated with a focus on the design of heuristic methods to solve different problem variants (Avella et al., 2004; Cornillier et al., 2008a, 2008b, 2009; Malépart et al., 2003; Ng et al., 2008; Taqa-allah et al., 2000).
3 Modelling Preliminaries
The demand structure, consumption rates, and penalty representations considered in this paper are similar to those used in Sherali and Al-Yakoob (2006a, 2006b), and are briefly stated below for the sake of completeness. Also, new notation related to the present work is introduced in this section.
3.1 Time-Horizon, Vessels, and Production Capacity
Let H be the set of days in the time-horizon, and let T be the set of types of vessels in the company’s fleet, where ${\Omega _{t}}$ represents the capacity of a vessel of type $t\in T$. For each $t\in T$, let $n=1,\dots ,{M_{t}}$ index the vessels of type t, and let ${O_{t}}$ and ${\mathit{CH}_{t}}={M_{t}}-{O_{t}}$ respectively, denote the number of company-owned vessels and the number of available vessels of this type that can be possibly chartered. The company-owned vessels are indexed by $n=1,\dots ,{O_{t}}$, while the chartered vessels are indexed by $n={O_{t}}+1,\dots ,{O_{t}}+{\mathit{CH}_{t}}\equiv {M_{t}}$. Accordingly, let ${\mathrm{\$ }_{t,n}}$ be the cost (in US dollars) of chartering a vessel n of type t, for each $n={O_{t}}+1,\dots ,{O_{t}}+{\mathit{CH}_{t}}$, and $t\in T$. Let ${\mathit{UT}_{t,n}}$ be the maximum number of days vessel n of type t can be used during the time-horizon, where this time restriction is typically needed for maintenance purposes. Let ${H^{t,n}}\subseteq \{1,\dots ,H\}$ be a subset of the time-horizon during which the ship n of type t will be available for use, where ${h_{1}^{t,n}}$ represents the first day in ${H_{}^{t,n}}$. A production capacity or a certain imposed supply quota for the product at the source is also specified, as given by Q.
3.2 Demand Structure and Penalty Functions
The demand structure of the problem is determined based on the following factors: (a) the storage capacity at the destination; (b) the initial level of each destination storage; (c) the rates of consumption at the destination; and (d) the production capacity at the source. The rates of consumption at a destination may vary due to seasonal considerations; however, we assume that this information is known a priori. For any given day of the time-horizon, there are minimum and maximum allowable limits imposed on the storage levels at the destinations, where certain agreed-upon penalties are incurred for violating these limits. At the beginning of the time-horizon, the storage level at the destination is given by w. This level may pertain to either a single storage unit or a collection of storage units at the destination; however, for the sake of modelling, we only deal with a combined aggregate storage level. Let SL${_{1}}$ and SL${_{2}}$ denote the minimum and maximum desired storage levels, respectively, at the client’s destination, which should be maintained to the extent possible in order to avoid penalties. Accordingly, let ${\Pi _{1}}>0$ and ${\Pi _{2}}>0$, respectively, denote the daily penalties for each shortage or excess unit at the destination, where ${\Pi _{1}}>{\Pi _{2}}$ is to emphasize the fact that shortage levels are more critical than excess levels. The permitted shortage and excess quantities at the destination with respect to the desired levels SL${_{1}}$ and SL${_{2}}$ are specified by the amounts ${A_{1}}$ and ${A_{2}}$, respectively. Let ${b_{1}}\equiv {\mathit{SL}_{1}}-{A_{1}}$ and ${b_{2}}\equiv {\mathit{SL}_{2}}+{A_{2}}$, and let $\mathit{UB}>{b_{2}}$ be a sufficiently large upper bound on the maximum allowable storage level on any given day of the time-horizon. Storage levels falling below ${b_{1}}$ or in excess of ${b_{2}}$ (up to $\mathit{UB}$), while permitted, are highly undesirable, and incur significantly greater penalties, respectively, given by ${\lambda _{1}}>{\Pi _{1}}$, and ${\lambda _{2}}>{\Pi _{2}}$ per unit.
Let ${R_{j}}$ denote the expected consumption rate at the destination on day j, for $j=1,\dots ,H$. The different daily consumption rates arise from possible seasonal changes during the time-horizon, as well as from client-specific considerations. Thus, the total cumulative consumption at the destination over the days $j=1,\dots ,h$ is given by ${\mathit{TC}_{h}}={\textstyle\sum _{j=1}^{h}}{R_{j}}$. The daily storage levels determine the overall penalty over the time-horizon, being given by the sum of Type I and Type II penalties as defined below, where ${S_{h}}$ denotes the storage level on day h and where ${\Pi _{1}}>0$, ${\Pi _{2}}>0$, ${\lambda _{1}}>{\Pi _{1}}$ and ${\lambda _{2}}>{\Pi _{2}}$ are as defined above.
Type I penalty: ${P_{I}}({S_{h}})={\Pi _{1}}$ maximum $\{0,({\mathit{SL}_{1}}-{S_{h}})\}$ if ${S_{h}}\in [{b_{1}},S{L_{1}}),$
$\hspace{54.06006pt}{P_{I}}({S_{h}})={\Pi _{2}}\hspace{2.5pt}\text{maximum}\hspace{2.5pt}\{0,({S_{h}}-S{L_{2}})\}\hspace{1em}\text{if}\hspace{2.5pt}{S_{h}}\in ({\mathit{SL}_{2}},{b_{2}}];\hspace{2.5pt}\text{and}$
Type II penalty: ${P_{\mathit{II}}}({S_{h}})=\left\{\begin{array}{l@{\hskip4.0pt}l}{\Pi _{1}}{A_{1}}+{\lambda _{1}}({b_{1}}-{S_{h}})\hspace{1em}& \text{if}\hspace{2.5pt}{S_{h}}\in [0,{b_{1}}),\\ {} {\Pi _{2}}{A_{2}}+{\lambda _{2}}({S_{h}}-{b_{2}})\hspace{1em}& \text{if}\hspace{2.5pt}{S_{h}}\in ({b_{2}},\mathit{UB}].\end{array}\right.$
Note that if ${S_{h}}\in [S{L_{1}},S{L_{2}}]$, then the storage level lies within the desired bounds and no penalty is induced. If ${S_{h}}\in [{b_{1}},{\mathit{SL}_{1}})\cup ({\mathit{SL}_{2}},{b_{2}}]$, then a Type I penalty is incurred based on the respective shortage or excess quantity. On the other hand, if ${S_{h}}\in [0,{b_{1}})\cup ({b_{2}},\mathit{UB}]$, then a relatively larger additional Type II penalty rate is imposed continuously beyond that of ${P_{I}}(.)$ to indicate the undesirabity of such a storage level on any given day of the time-horizon.
Proposition 1.
Let
where
Define the linear penalty function $P({S_{h}})={\Pi _{1}}{S_{2,h}}+{\Pi _{2}}{S_{4,h}}+{\lambda _{1}}{S_{3,h}}+{\lambda _{2}}{S_{5,h}}$, for $0\leqslant {S_{h}}\leqslant \mathit{UB}$. Then any minimization objective formulation that incorporates the term $P({S_{h}})$ defined above along with (1) and (2) will automatically enforce the sum of the Type I and Type II penalties ${P_{I}}({S_{h}})+{P_{\mathit{II}}}({S_{h}})$.
(2)
\[\begin{array}{l}\displaystyle {\mathit{SL}_{1}}\leqslant {S_{1,h}}\leqslant {\mathit{SL}_{2}},\hspace{1em}0\leqslant {S_{2,h}}\leqslant {A_{1}},\hspace{1em}0\leqslant {S_{3,h}}\leqslant {b_{1}},\\ {} \displaystyle 0\leqslant {S_{4,h}}\leqslant {A_{2}},\hspace{1em}\textit{and}\hspace{1em}0\leqslant {S_{5,h}}\leqslant UB-{b_{2}}.\end{array}\]Example.
The following simple example illustrates the Type I and Type II penalty structures. The various destinations storage levels associated with such penalties are illustrated in Fig. 1. Let $H=45$, where the per unit Type I and Type II penalties are given by $\{{\Pi _{1}}=100,{\Pi _{2}}=50\}$ and $\{{\lambda _{1}}=400,{\lambda _{2}}=150\}$, respectively. Assume that on day $h=1$ of the time horizon, the initial storage level at the destination is $w=5000$. Then the storage level on day h of the time horizon $({S_{h}})$ is given by 5000 plus the total amount delivered to the destination from day 1 to day h minus ${\mathit{TC}_{h}}$, where ${\mathit{TC}_{h}}$ denotes the total cumulative consumption at the destination over the days $j=1,\dots ,h$. Let ${R_{j}}=500$ for $j=1,\dots ,H$. Then ${\mathit{TC}_{h}}={\textstyle\sum _{j=1}^{h}}{R_{j}}=500h$, for $j=1,\dots ,h$.
Based on this information, we compute the incurred penalties for the days $h=5$, 15, 22, 29, and 36, assuming that the total amount delivered to the destination on days $h=5$, 15, 22, 29, and 36 are respectively given by 1000, 3000, 8000, 20000 and 25000. Note that, ${S_{5}}=3500$, ${S_{15}}=500$, ${S_{22}}=2000$, ${S_{29}}=10500$, and ${S_{36}}=12000$. The total Type I and Type II penalties incurred for each of the above five days, calculated based on Proposition 1, are given as follows:
\[\begin{array}{l}\displaystyle {P_{I}}({S_{5}})={P_{\mathit{II}}}({S_{5}})=0\hspace{2.5pt}(\text{no Penalty}),\hspace{2em}{P_{I}}({S_{15}})=0,\hspace{2em}{P_{\mathit{II}}}({S_{15}})=400000,\\ {} \displaystyle {P_{I}}({S_{22}})=25000,\hspace{2em}{P_{\mathit{II}}}({S_{22}})=0,\hspace{2em}{P_{I}}({S_{29}})=25000,\hspace{2em}{P_{\mathit{II}}}({S_{29}})=0,\\ {} \displaystyle {P_{I}}({S_{36}})=0,\hspace{1em}\text{and}\hspace{1em}{P_{\mathit{II}}}({S_{36}})=200000.\end{array}\]
3.3 Viable Trips, Average Speeds, and Transportation Costs
There are five permissible journey itineraries, denoted by ${J_{j}}$, for $j=1,\dots ,5$, as defined next. Let $i=1,2,3$ index the source, transshipment facility, and destination, respectively, and let $\mathit{SFD}=\{1,2,3\}$. Then the permissible journeys are given as follows:
${J_{1}}=(1,3,1)$: source → destination → source;
${J_{2}}=(1,2,1)$: source → transshipment facility → source;
${J_{3}}=(1,3,2)$: source → destination → transshipment facility;
${J_{4}}=(2,3,1)$: transshipment facility → destination → source;
${J_{5}}=(2,3,2)$: transshipment facility → destination → transshipment facility.
Let $J={\textstyle\bigcup _{j=1}^{5}}{J_{j}}$, which is composed of all permissible journeys. For $j=1,\dots ,5$, let ${J_{j}}=({p_{1}^{j}},{p_{2}^{j}},{p_{3}^{j}})$ characterize the itinerary of Journey ${J_{j}}$ as identified above, where ${p_{1}^{j}},{p_{2}^{j}}$, and ${p_{3}^{j}}$ respectively represent the loading, unloading, and termination points of Journey ${J_{j}}$, and let ${J_{j}^{1}}\equiv ({p_{1}^{j}},{p_{2}^{j}})$ and ${J_{j}^{2}}\equiv ({p_{2}^{j}},{p_{3}^{j}})$ denote the two sequential legs of this itinerary. For $i=1,2$, let ${I_{i}^{1}}$ be a subset of J such that i is the loading point of any journey in ${I_{i}^{1}}$. Hence, ${I_{1}^{1}}=\{{J_{1}},{J_{2}},{J_{3}}\}$ and ${I_{2}^{1}}=\{{J_{4}},{J_{5}}\}$. For $i=2,3$, let ${I_{i}^{2}}$ be a subset of J such that i is the intermediate point (unloading point) of any journey in ${I_{i}^{2}}$. Thus, ${I_{2}^{2}}=\{{J_{2}}\}$ and ${I_{3}^{2}}=\{{J_{1}},{J_{3}},{J_{4}},{J_{5}}\}$. For $i=1,2$, let ${I_{i}^{3}}$ be a subset of J such that i is the terminal point of any journey in ${I_{i}^{3}}$. Accordingly, ${I_{1}^{3}}=\{{J_{1}},{J_{2}},{J_{4}}\}$ and ${I_{2}^{3}}=\{{J_{3}},{J_{5}}\}$.
Remark 1.
We assume that partial loading/unloading is not allowed and that all journeys commence and terminate at either the source or the transshipment facility. Partial loading/unloading could, however, further enhance the overall operational efficiency and cost by allowing other types of journeys in the operation, and could be likewise handled in a similar fashion. We defer this for future research.
The average speed (knots per hour) of a vessel of type t is given by ${V_{t}^{F}}$ when the vessel is fully loaded, and ${V_{t}^{E}}$ when the vessel is empty. The permitted average daily utilization of a vessel of type t (in hours per day) is given by ${\mathit{DU}_{t}^{F}}$ when the vessel is fully loaded, and ${\mathit{DU}_{t}^{E}}$ when the vessel is empty. Let ${C_{t}^{\mathit{PU},F}}$and ${C_{t}^{\mathit{PU},E}}$ be the daily transportation costs associated with a vessel of type t when the vessel is fully loaded and empty, respectively. Note that ${C_{t}^{\mathit{PU},F}}$ and ${C_{t}^{\mathit{PU},E}}$ are calculated as the product of the vessel capacity ${\Omega _{t}}$ and the per unit transportation cost when the vessel is fully loaded and empty, respectively.
Next, we compute the cost associated with journeys in the set J. For let, denote the distance (in knots) between α and β. For $j\in \{1,\dots ,5\}$, let $\mathit{TD}({J_{j}})$ be the total distance that a vessel travels to complete Journey ${J_{j}}$. Accordingly, we have
\[\begin{array}{l}\displaystyle \mathit{TD}({J_{1}})=\delta (1,3)+\delta (3,1),\hspace{2em}\mathit{TD}({J_{2}})=\delta (1,2)+\delta (2,1),\\ {} \displaystyle \mathit{TD}({J_{3}})=\delta (1,3)+\delta (3,2),\hspace{2em}\mathit{TD}({J_{4}})=\delta (2,3)+\delta (3,1),\hspace{1em}\text{and}\\ {} \displaystyle \mathit{TD}({J_{5}})=\delta (2,3)+\delta (3,2).\end{array}\]
Let ${\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}=\frac{\delta ({p_{1}^{j}},{p_{2}^{j}})}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}}$ and ${\mathit{ND}_{t}^{({p_{2}^{j}},{p_{3}^{j}})}}=\frac{\delta ({p_{2}^{j}},{p_{3}^{j}})}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}}$, where the former represents the average travel time (in days) required from point ${p_{1}^{j}}$ to point ${p_{2}^{j}}$, and the latter represents the average travel time (in days) required from point ${p_{2}^{j}}$ to point ${p_{3}^{j}}$. Also, let ${\mathit{ND}_{t}^{j}}={\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}+{\mathit{ND}_{t}^{({p_{2}^{j}},{p_{3}^{j}})}}$.Hence, the transportation cost of a vessel of type t that undertakes journey ${J_{j}}$ is given by ${C_{t,j}}={C_{t,j}^{F}}+{C_{t,j}^{E}}$, where ${C_{t,j}^{F}}={C_{t}^{\mathit{PU},F}}{\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}$ and ${C_{t,j}^{E}}={C_{t}^{\mathit{PU},E}}{\mathit{ND}_{t}^{({p_{2}^{j}},{p_{3}^{j}})}}$ represent the total transportation cost for a vessel of type t when this vessel is fully loaded and empty, respectively.
It is worth mentioning that the transportation costs of vessels were assumed to be fixed and known a priori in Sherali and Al-Yakoob (2006a, 2006b). In contrast, in this paper, the transportation costs are fixed for the journeys that do not involve the transshipment facility; otherwise, these costs are functions of the location of the transshipment facility.
3.4 Transshipment Facility
Let ${H^{f}}\subseteq \{1,\dots ,H\}$ be a contiguous subset of the time-horizon during which the transshipment facility will be available for leasing. Let ${h_{1}^{f}}$ be the first day in ${H^{f}}$ and let ${w^{f}}$ denote the storage level at the transshipment facility on Day ${h_{1}^{f}}$. Let $\mathit{CD}$ denote the capacity of the transshipment facility and let $\mathrm{\$ }\mathrm{\$ }$ be a fixed cost (in U.S. dollars) of leasing the transshipment facility during ${H^{f}}$, where we set $\mathrm{\$ }\mathrm{\$ }\equiv 0$ if the transshipment facility is owned by the company. Let MC be the daily maintenance cost of the transshipment facility, and so, $|{H^{f}}|$MC represents the total maintenance cost over the entire leasing duration. Let AR define the allowable region in which the transshipment facility may be located.
4 An Aggregated Mathematical Programming Model
Following Sherali and Al-Yakoob (2006a, 2006b), we initially formulated an exact disaggregated model (EM) for the problem under consideration. This model identified individual vessel characteristics related to capacity, travel time, maintenance requirements, etc. Our preliminary computational experiments revealed that Model EM was intractable for large test instances because of its relative size and due to the symmetry inherent within this model, in the sense that the vessel schedules of identically characterized vessels (same type, availability, and utilization requirements) could be permuted to generate reflections of essentially the same solution (related symmetry-defeating constraints Sherali and Smith (2001) were unhelpful in this regard). Hence, we propose an aggregated formulation (AF), which retains the salient features of EM while attempting to determine the aggregate number of vessels of each type to be dispatched on a prescribed trip from the source to the transshipment facility on a given day, rather than prescribing schedules of individual vessels, as done in EM. This aggregate representation results in a fewer number of variables and constraints than those in EM, and moreover, Model AF automatically suppresses the aforementioned symmetry effects that occur within Model EM. Hence, for the sake of brevity in presentation, we focus in this paper on Model AF, and combine this with the methods described in Sherali and Al-Yakoob (2006b) and Al-Yakoob and Sherali (2013) to design a hybrid, two-stage solution approach, as subsequently discussed in Section 6.
4.1 Variable Definition and Related Issues
Let ${\eta _{{h_{1}},t,j,h}^{1}}$ and ${\eta _{{h_{1}},t,j,h}^{2}}$ for $\forall {h_{1}},t,j,h>{h_{1}}$ be binary indicator parameters that take on a value of 1 if and only if the delivery date h of vessels of type t traverses journey ${J_{j}}$ starting on day h; exceeds the corresponding start time plus the journey duration, and serves to set the associated duration-incompatible u-variables that are defined subsequently to zero. Define ${x_{h,t,j}}$ to be an integer variable that represents the number of vessels of type t that traverse journey ${J_{j}}$ starting on day h; and let ${u_{{h_{1}},t,j,h}^{1}}={x_{{h_{1}},t,j}}{\eta _{{h_{1}},t,j,h}^{1}}$ and ${u_{{h_{1}},t,j,h}^{2}}={x_{{h_{1}},t,j}}{\eta _{{h_{1}},t,j,h}^{2}}$, $\forall {h_{1}},t,j,h>{h_{1}}$. Hence, ${u_{{h_{1}},t,j,h}^{1}}$ is an integer variable that represents the number of vessels of type t-th at start journey ${J_{j}^{}}$ on day ${h_{1}}$ and that arrive at point ${p_{2}^{j}}$ on or before day h. Similarly, ${u_{{h_{1}},t,j,h}^{2}}$ is an integer variable that represents the number of vessels of type t that start journey ${J_{j}}$ on day ${h_{1}}$ and that terminate this journey at point ${p_{3}^{j}}$ on or before day h.
Define ${y_{h,t,s}}$ to be an integer decision variable that represents the maximum number of vessels of type t that are available for dispatching from the source (if $s=1$) or from the transshipment facility (if $s=2$) on day h.
As mentioned earlier, vessels might become available for use at different days of the time-horizon due to, for example, their involvement in trips from previous demand contracts that will terminate sometime during the current time-horizon. Hence, we let ${O_{h,t}}$ be the number of self-owned vessels of type t that will become available for use for the first time at the source on day h of the time-horizon, and we let ${\mathit{CH}_{h,t}}$ be the number of vessels of type t that will become available for the first time at the source for chartering on day h of the time-horizon. Let ${\alpha _{h,t}}={O_{h,t}}+{\mathit{CH}_{h,t}}$ and note that ${O_{t}}={\textstyle\sum _{h}}{O_{h,t}}$ and ${\mathit{CH}_{t}}={\textstyle\sum _{h}^{}}{\mathit{CH}_{h,t}}$. Accordingly, we let ${z_{h,t}}$ be an integer variable that denotes the number of vessels of type t that are actually selected (from among the ${\mathit{CH}_{h,t}}$ vessels) for chartering on day h of the time-horizon, and we let ${\mathrm{\$ }^{h,t}}$ denote the average chartering cost of a vessel of type t that will become available for potential use on day h of the time-horizon.
Remark 2.
Tentatively, consider another binary variable W defined to take the value of one if the transshipment depot is leased during ${H^{f}}$, and is zero otherwise. Note that $W=0$ implies that the transshipment facility is not selected in the operation and hence only journey ${J_{1}}$ is permitted, whence the problem simplifies to that studied in Sherali and Al-Yakoob (2006b). Accordingly, for the sake of ease in modelling, we can consider two cases separately, namely, when $W=0$ and $W=1$, and compare the respective operational costs. Since the case $W=0$ has been studied in Sherali and Al-Yakoob (2006b), we focus in the remainder of this paper on the case when the transshipment facility (to-be-located) is leased $(W=1)$.
Let ${A_{h,t}}$ be the subset of indices for vessels of type t (both self-owned and vessels available for chartering) that will become available for use at the source for the first time on day h of the time-horizon. Hence, for a given day h and vessel type t, we let ${\mathit{UT}_{h,t}}=\frac{{\textstyle\sum _{n\in {A_{h,t}}}^{}}{\mathit{UT}_{t,n}}}{{\alpha _{h,t}}}$ represent the average usage allowance for a vessel of type t that will become available for use for the first time on day h of the time-horizon. Accordingly, ${\mathit{UT}_{t}}=\frac{{\textstyle\sum _{h}^{}}{\alpha _{h,t}}U{T_{h,t}}}{{\textstyle\sum _{h}^{}}{\alpha _{h,t}}}$ gives the average usage allowance for a vessel of type t.
Note that in this aggregate model, the variable ${y_{h,t,s}}$ represents the maximum number of vessels of type t that could be consigned on day h as necessary; the actual number used, and in particular the chartering decisions, are governed in this model via the dispatching variables ${x_{h,t,j}}$.
Let ${\phi _{x}}$, ${\phi _{{\eta ^{1}}}}$, ${\phi _{{\eta ^{2}}}}$, ${\phi _{{u^{1}}}}$, ${\phi _{{u^{2}}}}$, and ${\phi _{y}}$ be the set of indices of the x-, ${\eta ^{1}}$-, ${\eta ^{2}}$-, ${u^{1}}$-, ${u^{2}}$- and y-variables, respectively, that are a priori restricted to be zero, or fixed at some known positive integer values.
4.2 Model Constraints
The constraints of Model AF are formulated next below.
A) Examining arrival dates of shipments. Consider a journey ${J_{j}}\in J$ and a vessel of type $t\in T$. The following constraints assign proper binary values for the ${u^{1}}$-variables based on whether this vessel starts journey j on day ${h_{1}}$ and delivers its shipload at location ${p_{2}^{j}}$ on or before day h. Hence, Constraint (3) enforces the defined value of ${\eta _{{h_{1}},t,j,h}^{1}}$ to equal to 1 or 0 based on whether ${h_{1}}+{\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}-h\leqslant 0$ or ${h_{1}}+{\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}-h\geqslant 1$, respectively, and then Constraints (4)–(6) along with ${u_{{h_{1}},t,j,h}^{1}}\geqslant 0$ enforce the product relationship ${u_{{h_{1}},t,j,h}^{1}}={x_{{h_{1}},t,j}}{\eta _{{h_{1}},t,j,h}^{1}}$, $\forall {h_{1}},t,j,h>{h_{1}}$.
(3)
\[\begin{array}{l}\displaystyle 1-{\eta _{{h_{1}},t,j,h}^{1}}(h+1)\leqslant {h_{1}}+{\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}-h\leqslant \big(1-{\eta _{{h_{1}},t,j,h}^{1}}\big)(H-h),\\ {} \displaystyle \hspace{1em}\forall {h_{1}},t,j,h>{h_{1}},\end{array}\](4)
\[ {u_{{h_{1}},t,j,h}^{1}}\leqslant {x_{{h_{1}},t,j}},\hspace{1em}\forall {h_{1}},t,j,h>{h_{1}},\]
B) Examining arrival dates of vessels. Consider a journey ${J_{j}}\in J$ and a vessel of type t. Similar to the above constraints for the ${u^{1}}$-variable, the following constraints assign proper integer values for the ${u^{2}}$-variables based on whether this vessel starts journey j on day ${h_{1}}$ and returns to location ${p_{3}^{j}}$ on or before day h, where recall that ${u_{{h_{1}},t,j,h}^{2}}={x_{{h_{1}},t,j}}{\eta _{{h_{1}},t,j,h}^{2}}$, $\forall {h_{1}},t,j,h>{h_{1}}$. Accordingly, Constraint (7) enforces the defined binary value of ${\eta _{{h_{1}},t,j,h}^{2}}$ and then (8)–(10) (along with ${u_{{h_{1}},t,j,h}^{2}}\geqslant 0)$ represent the stated product relationship.
(7)
\[\begin{array}{l}\displaystyle 1-{\eta _{{h_{1}},t,j,h}^{2}}(h+1)\leqslant {h_{1}}+{\mathit{ND}_{t}^{j}}-h\leqslant \big(1-{\eta _{{h_{1}},t,j,h}^{2}}\big)(H-h),\\ {} \displaystyle \hspace{1em}\forall {h_{1}},t,j,h>{h_{1}},\end{array}\]
C) Representation of the destination’s storage level and related penalties. The daily storage level of the product at the destination must remain within $[{\mathit{SL}_{1}},{\mathit{SL}_{2}}]$ to the extent possible, where appropriate daily penalties are imposed otherwise based on the specific levels of the storage as discussed above (see Proposition 1). The representation of the storage level of the destination is given by the following constraints:
Constraint (11) gives the storage level on day h based on the daily consumption rates and the shipments of the product that are delivered by day h. Constraint (12) represents ${S_{h}}$ in terms of ${S_{1,h}},{S_{2,h}},{S_{3,h}},{S_{4,h}}$, and ${S_{5,h}}$ as discussed above in Proposition 1, where (2) holds, $\forall h$, so that appropriate penalties would be incurred in the objective function based on this representation.
D) Representation of the transshipment facility’s storage level. The following constraint provides a representation of the daily storage level ${L_{h}}$ for day h at the transshipment facility, where ${\mathit{LB}^{f}}$ and ${\mathit{UB}^{f}}$ are some specified lower and upper bounds on the storage level of the facility:
along with ${\mathit{LB}^{f}}\leqslant {L_{h}}\leqslant {\mathit{UB}^{f}}$, $\forall h\in {H^{f}}$.
(13)
\[ {L_{h}}={w^{f}}+\sum \limits_{{h_{1}^{f}}\leqslant {h_{1}}<h}\sum \limits_{t}\sum \limits_{j\in {I_{2}^{2}}}{\Omega _{t}}{u_{{h_{1}},t,j,h}^{1}}-\sum \limits_{{h_{1}^{f}}\leqslant {h_{1}}\leqslant h}\sum \limits_{t}\sum \limits_{j\in {I_{2}^{1}}}{\Omega _{t}}{x_{{h_{1}},t,j}},\hspace{1em}\forall h\in {H^{f}},\]The storage level of the transshipment facility on a day $h\in {H^{f}}$ is determined based on the initial level ${w^{f}}$ of the facility on day ${h_{1}^{f}}$, plus the total amount of the product that is delivered to the facility by day h, and minus the total amount of the product that is shipped from the facility by day h.
E) Availability of vessels and vessel chartering. The vessel availability constraints are given as follows:
(14)
\[\begin{array}{l}\displaystyle {y_{h,t,s}}={y_{h-1,t,s}}-\sum \limits_{j\in {I_{s}^{1}}}{x_{h-1,t,j}}+\sum \limits_{{h_{1}}<h}\sum \limits_{j\in {I_{s}^{1}}}\big({u_{{h_{1}},t,j,h}^{2}}-{u_{{h_{1}},t,j,(h-1)}^{2}}\big)+{O_{h,t}}+{z_{h,t}},\\ {} \displaystyle \hspace{1em}h\geqslant 2,\hspace{2.5pt}t,s=1,2,\end{array}\]A vessel of type t can be dispatched on a journey beginning from point $s\in \{1,2\}$ on day h only if it is available at s on that day. This vessel is available at s on day h if either the vessel was available at s on the previous day and it was not dispatched, or this vessel was not available there during the previous day, but it became available on the current day. On the other hand, this vessel is unavailable on day h at location s if it was available there on the previous day and it was dispatched on that day (noting that ${N_{t}^{j}}\geqslant 2$ for all vessel-types t), or it was unavailable on the previous day and it did not arrive on the current day. Constraint (14) examines the availability of vessel s of type t at location s on day h by considering these cases, and then Constraint (15) permits the dispatching of vessels conditioned upon their availability. Note that the hand-side of Constraint (14) accounts for the first-time availabilities of the self-owned and chartered vessels via the last two terms, where Constraint (16) serves to initialize the y-variables for $h=1$ (noting that ${y_{1,t,2}}=0$, $\forall t$). Furthermore, (17) is a definitional constraint that computes the total number of vessels of type t that are actually selected for chartering.
F) Capacity restrictions and maintenance requirements. A production capacity or certain imposed daily supply quota Q restricts the maximum cumulative amount of the product that can be shipped from the source to the transshipment facility and to the destination over the periods $1\leqslant {h_{1}}\leqslant h$, $\forall h$. This restriction is represented by the following constraint:
Furthermore, for maintenance purposes, the following constraint enforces that vessel of type t can be used for at most ${\mathit{UT}_{t}}$ days during the time-horizon:
Other forms of scheduled maintenance restrictions can be also accommodated in the model by setting certain x-variables to zero. Note that if $j\in \{2,3,4,5\}$, i.e. ${J_{j}}$ is a journey that involves the transshipment facility, then $[{\mathit{ND}_{t}^{j}}]{x_{h,t,n,j}}$ is a product of two variables.
(18)
\[ \sum \limits_{{h_{1}}\leqslant h}\sum \limits_{t}\sum \limits_{j\in {I_{1}^{1}}}{\Omega _{t}}{x_{{h_{1}},t,j}}\leqslant hQ,\hspace{1em}\forall h.\](19)
\[ \sum \limits_{h}\sum \limits_{j}\big[{\mathit{ND}_{t}^{j}}\big]{x_{h,t,j}}\leqslant {\mathit{UT}_{t}}({O_{t}}+{z_{t}}),\hspace{1em}\forall t.\]4.3 Objective Function and Overall Model
The objective function is composed of the following terms:
a) Operational costs of vessels (both self-owned and chartered):
b) Penalties resulting from violating the allowable storage level $[{\mathit{SL}_{1}},{\mathit{SL}_{2}}]$ at the destination, which are given as follows based on Proposition 1:
c) Chartering expenses of vessels: ${\textstyle\sum _{h}}{\textstyle\sum _{t}}$${\mathrm{\$ }^{h,t}}{z_{h,t}}$.
d) Leasing and maintenance expenses of the transshipment facility: based on Remark 2, a fixed cost given by $\mathrm{\$ }\mathrm{\$ }+\| {H^{f}}\| \mathit{MC}$ is incurred in the objective function.
The objective function terms (a)–(d) along with the constraints formulated above yield the following model AF, where all indices are assumed to take on only their respective relevant values:
\[\begin{array}{l}\displaystyle \textbf{AF}:\hspace{2.5pt}\text{Minimize}\hspace{2.5pt}\sum \limits_{h}{\sum \limits_{t}^{}}\sum \limits_{j}{C_{t,j}}{x_{h,t,j}}+\sum \limits_{h}{\Pi _{1}}{S_{2,h}}+{\sum \limits_{h}^{}}{\Pi _{2}}{S_{4,h}}\\ {} \displaystyle \hspace{1em}+{\sum \limits_{h}^{}}{\lambda _{1}}{S_{3,h}}+{\sum \limits_{h}^{}}{\lambda _{2}}{S_{4,h}}+\sum \limits_{h}\sum \limits_{t}{\mathrm{\$ }^{h,t}}{z_{h,t}},\end{array}\]
subject to (3)–(19).Where,
${x_{h,t,j}}\in \{0,1,\dots ,{M_{t}}\}$, $\forall (h,t,j)\notin {\phi _{x}}$, and fixed at zero otherwise.
${\eta _{{h_{1}},t,j,h}^{1}}\in \{0,1\}$, $\forall ({h_{1}},t,j,h>{h_{1}})\notin {\phi _{{\eta ^{1}}}}$, and fixed at zero or one otherwise.
${\eta _{{h_{1}},t,j,h}^{2}}\in \{0,1\}$, $\forall ({h_{1}},t,j,h>{h_{1}})\notin {\phi _{{\eta ^{2}}}}$, and fixed at zero or one otherwise.
${u_{{h_{1}},t,j,h}^{1}}\geqslant 0$, $\forall ({h_{1}},t,j,h>{h_{1}})\notin {\phi _{{u^{2}}}}$, and fixed at zero otherwise.
${u_{{h_{1}},t,j,h}^{2}}\geqslant 0$, $\forall ({h_{1}},t,j,h>{h_{1}})\notin {\phi _{{u^{2}}}}$, and fixed at zero otherwise.
${y_{h,t,s}}\in [0,{M_{t}}]$, $\forall (h\geqslant 2$, $t,s\in \{1,2\})\notin {\phi _{y}}$, and fixed at zero or a positive integer otherwise.
${z_{h,t}}\in \{1,\dots ,{\mathit{CH}_{h,t}}\}$, $\forall h,t$.
$S{L_{1}}\leqslant {S_{1,h}}\leqslant S{L_{2}}$, $0\leqslant {S_{2,h}}\leqslant {A_{1}}$, $0\leqslant {S_{3,h}}\leqslant {b_{1}}$, $0\leqslant {S_{4,h}}\leqslant {A_{2}}$.
$0\leqslant {S_{5,h}}\leqslant \mathit{UB}-{b_{2}}$, $\forall h$.
$L{B^{f}}\leqslant {L_{h}}\leqslant {\mathit{UB}^{f}}$, $\forall h\in {H^{f}}$.
Noting that the transshipment facility location must belong to AR.
Proposition 2.
If we enforce binary restrictions on the $(x,{\eta ^{1}},{\eta ^{2}})$-variables, then the $({u_{}^{1}},{u_{}^{2}},y,z)$-variables will automatically turn out to be binary-valued when simply restricted to be continuous.
Proof.
Suppose that x, ${\eta ^{1}}$, and ${\eta ^{2}}$ are binary-valued, and that ${u^{1}}$, ${u^{2}}$, y, and z are restricted to be continuous as specified in Model AF. Then the $({u^{1}},{u^{2}})$-variables will be binary-valued by Constraints (4)–(6) and Constraints (8)–(10), respectively. This implies that the y-variables will also automatically be binary-valued by Constraint (14), noting that all the y-variables corresponding to $h=1$ are initialized at zero or one. This follows from the recursive relation in (14) that defines ${y_{h,t,s}}$ for $h\geqslant 2$ in terms of ${y_{h-1,t,s}}$ and a subset of the x- and u-variables, whereby the y-variables will be integral within $[0,1]$. Moreover, the binary restrictions on the z- variables will then automatically hold via Constraint (16), noting that these variables have a positive coefficient in the objective function. □
5 Transshipment Facility Location and Related Costs
This section discusses modelling considerations related to the transshipment facility, which are instrumental in ensuring the tractability of Model AF. In this context, we examine some plausible line segments along shore lines of admissible land bases where a transshipment facility could be possibly located. Note that this could also include some off-shore sites, whenever the location of a storage facility rig is feasible at such locations. Hence, we consider a suitable collection of line segments and solve Model AF for each such segment in order to determine a set of cost-effective locations for the storage facility, one for each line segment, and we then use the column generation method of Al-Yakoob and Sherali (2013) with this set of discretized locations to generate improved cost-effective vessel schedules. This hybrid, two stage algorithm is proposed in the next section.
Accordingly, let ${l_{1}},\dots ,{l_{N}}$ represent the set of such (directed) plausible line segments. For $i=1,\dots ,N$, let AF$({l_{i}})$ denote Model AF with ${l_{i}}$ expressing the set of feasible potential locations for the transshipment facility. Letting $v(\mathrm{P})$ denote the optimal objective function value of any model P, we have that $v(\mathrm{AF})={\min _{i=1\dots ,N}}\{v(\mathrm{AF}({l_{i}})\}$ gives the optimal value for Model AF. Let l be any of the directed segments from $\{{l_{1}},\dots ,{l_{N}}\}$ and let ${l^{0}}$ and ${l^{1}}$ respectively denote the initial and final points of l. For any $\xi \in [0,1]$, let $p(\xi )$ be the point on l that represents a convex combination of ${l^{0}}$ and ${l^{1}}$, with respectively weights ξ and $(1-\xi )$. Hence, the distance from the source to $p(\xi )$ is given by $\delta (1,p(\xi ))=\delta (1,{l^{0}})+\delta ({l^{0}},p(\xi ))$, where $\delta (1,{l^{0}})$ is known a priori and $\delta ({l^{0}},p(\xi ))=(1-\xi )\delta ({l^{0}},{l^{1}})$. Also, the distance from $p(\xi )$ to the destination is given by $\delta (p(\xi ),3)=\delta (p(\xi ),{l^{1}})+\delta ({l^{1}},3)$, where $\delta (p(\xi ),{l^{1}})=\xi \delta ({l^{0}},{l^{1}})$ and $\delta ({l^{1}},3)$ is known a priori.
Recall that
\[ {\mathit{ND}_{t}^{({p_{1}^{j}},{p_{2}^{j}})}}=\frac{\delta ({p_{1}^{j}},{p_{2}^{j}})}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}},\]
which is fixed for $j\in \{1,3\}$, but is a function of the location of the transshipment facility for $j\in \{2,4,5\}$. Likewise,
\[ {\mathit{ND}_{t}^{({p_{2}^{j}},{p_{3}^{j}})}}=\frac{\delta ({p_{2}^{j}},{p_{3}^{j}})}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}}\]
is fixed for $j\in \{1,4\}$, and is a function of the location of the transshipment facility for $j\in \{2,3,5\}$. Hence, given that the transshipment facility is located at $p(\xi )$ for some $\xi \in [0,1]$, we have that
\[\begin{array}{l}\displaystyle {\mathit{ND}_{t}^{(1,2)}}=\frac{\delta (1,2)}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}}=\frac{\delta (1,{l^{0}})+(1-\xi )\delta ({l^{0}},{l^{1}})}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}},\\ {} \displaystyle {\mathit{ND}_{t}^{(2,1)}}=\frac{\delta (1,2)}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}}=\frac{\delta (1,{l^{0}})+(1-\xi )\delta ({l^{0}},{l^{1}})}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}},\\ {} \displaystyle {\mathit{ND}_{t}^{(2,3)}}=\frac{\delta (2,3)}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}}=\frac{\xi \delta ({l^{0}},{l^{1}})+\delta ({l^{1}},3)}{{V_{t}^{F}}{\mathit{DU}_{t}^{F}}},\\ {} \displaystyle {\mathit{ND}_{t}^{(3,2)}}=\frac{\delta (2,3)}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}}=\frac{\xi \delta ({l^{0}},{l^{1}})+\delta ({l^{1}},3)}{{V_{t}^{E}}{\mathit{DU}_{t}^{E}}}.\end{array}\]
Hence, ${C_{t,2}^{F}}={C_{t}^{\mathit{PU},F}}{\mathit{ND}_{t}^{(1,2)}}$, ${C_{t,4}^{F}}={C_{t,5}^{F}}={C_{t}^{\mathit{PU},F}}{\mathit{ND}_{t}^{(2,3)}}$, ${C_{t,2}^{E}}={C_{t}^{\mathit{PU},E}}{\mathit{ND}_{t}^{(2,1)}}$, and ${C_{t,3}^{E}}={C_{t,5}^{E}}={C_{t}^{\mathit{PU},E}}={\mathit{ND}_{t}^{(3,2)}}$.Thus the total operational cost of vessels in Model AF can be written as follows:
\[\begin{aligned}{}\sum \limits_{h}\sum \limits_{t}\sum \limits_{j}{C_{t,j}}{x_{h,t,j}}=& \sum \limits_{h}\sum \limits_{t}\sum \limits_{j}{C_{t,j}^{F}}{x_{h,t,j}}+\sum \limits_{h}\sum \limits_{t}\sum \limits_{j}{C_{t,j}^{E}}{x_{h,t,j}}\\ {} =& \sum \limits_{h}\sum \limits_{t}\big[{C_{t,1}^{F}}{x_{h,t,1}}+{C_{t,3}^{F}}{x_{h,t,3}}\big]\\ {} & +\sum \limits_{h}\sum \limits_{t}\big[{C_{t,1}^{E}}{x_{h,t,1}}+{C_{t,4}^{E}}{x_{h,t,4}}\big]\\ {} & +\sum \limits_{h}\sum \limits_{t}\big[{C_{t,2}^{F}}{x_{h,t,2}}+{C_{t,4}^{F}}{x_{h,t,4}}+{C_{t,5}^{F}}{x_{h,t,5}}\big]\\ {} & +\sum \limits_{h}\sum \limits_{t}\big[{C_{t,2}^{E}}{x_{h,t,2}}+{C_{t,3}^{E}}{x_{h,t,4}}+{C_{t,5}^{E}}{x_{h,t,5}}\big].\end{aligned}\]
The product relationships of the type $\xi x$ in AF, can be linearized as described next, where we drop the indices on the x-variables for ease in presentation.
For Model AF, we define $\hat{x}=\xi x$, and accordingly enforce this product relationship using the special structured Reformulation-Linearization Technique of Sherali et al. (1998) as follows:
where
where
Note that Constraints (20)–(21) represent the integral x-variables in terms of the binary λ-variables, and then Constraints (23)–(24) enforce the required product relationship $\hat{x}=\xi x$, noting that ${\hat{x}_{k}}\equiv \xi {\lambda _{k}}$, for $k=0,1,\dots ,M.$ The nonlinear terms that appear in Constraints (19) of AF for the case of $j\in \{2,3,4,5\}$ can be linearized similarly as above.
6 Algorithmic Approach
In this section, we design a two-stage hybrid algorithm that combines the approaches in Sherali and Al-Yakoob (2006b) and Al-Yakoob and Sherali (2013) along with Model AF of this paper to generate good quality solutions for our underlying problem. Note that in Sherali and Al-Yakoob (2006b), we formulated exact and aggregated vessel scheduling models (denoted by VSM and AVSM, respectively) for the case of multiple sources and multiple destinations, while considering the leasing of transshipment storage facilities, each of which has a known fixed location. The model AVSM decides on the number of vessels of each type to be dispatched on days of the time-horizon while Model VSM determines the dispatching of individual vessels. Model AVSM can be solved more efficiently than Model VSM because the latter involves a relatively large number of variables and also suffers from symmetry effects (see Sherali and Smith, 2001). For practical sized test problems that could not be solved directly via Model AVSM, we proposed in Sherali and Al-Yakoob (2006b) a rolling horizon algorithm, denoted by RHA, which generates solutions for Model AVSM in an iterative fashion by partitioning the time-horizon, where at each iteration, the integrality of only a subset of the discrete variables pertaining to imminent time periods is enforced, while certain other discrete variables are fixed as determined from previous iterations, and where the remaining variables are relaxed to take on continuous values in their respective domains. The values of the integer-restricted variables are then fixed as obtained from the resulting solution, and the horizon is accordingly rolled forward for the next iteration. In Al-Yakoob and Sherali (2013), a specially designed MIP model (VM) was proposed for the problem studied in Sherali and Al-Yakoob (2006b), while allowing the transshipment facility to be situated at a point selected from a discretized ad hoc list of potential locations. A column generation heuristic (CGH) was then devised in Al-Yakoob and Sherali (2013) to solve Model VM. Computational results revealed that CGH substantially outperformed the approach of using RHA to solve AVSM with respect to both the quality of solution obtained and the CPU time consumed, even when each transshipment facility was assigned a single fixed location in the former approach. Moreover, the results in Al-Yakoob and Sherali (2013) demonstrated that the consideration of multiple alternative locations for each available transshipment facility is likely to further reduce costs.
Based on the foregoing observations, we propose next a two-stage hybrid algorithm, denoted by HA, for our underlying problem. In the first stage of HA, we apply Algorithm RHA of Sherali and Al-Yakoob (2006b) to solve Model AF based on a specified set of plausible line segments, as discussed in Section 5, in order to determine a discretized set of best-cost potential facility locations, up to one for each line segment. Let this set be denoted by BL. In determining BL, note that for any given line segment, the solution of AF may determine a location for the transshipment facility if its use at this selected location reduces the overall operational cost. In the second stage of HA, instead of simply selecting the best location within BL by examining the objective values of these locations as given via Model AF, we formulate Model VM and solve it using CGH as in Al-Yakoob and Sherali (2013), while allowing the facility to be situated at any location from the set BL.
7 Computational Results and Related Issues
We consider 10 test problem instances based on simulated data pertaining to the Kuwait Petroleum Corporation (KPC), denoted by ${P_{1}},\dots ,{P_{10}}$. Detailed data for these test problems are available at www.al-yakoob.com under “Research Related Statistics”. Aside from using different combinations of time horizons, owned and chartered vessels, and travel times, the basic data for these test problems is essentially the same as that of Sherali and Al-Yakoob (2006a).
The following notation is used in presenting our computational results:
-
• $\overline{\mathrm{M}}$: Linear Programming (LP) relaxation of model M.
-
• RHA: The rolling horizon algorithm described in Sherali and Al-Yakoob (2006b), which is used in the present paper for solving Model AF.
-
• VM: A specially designed MIP model proposed for the problem studied in Sherali and Al-Yakoob (2006b), while allowing the transshipment facility to be situated at a point selected from a discretized ad hoc list of potential locations.
-
• HA: The two-stage hybrid algorithm described in Section 6.
-
• CGM: The Column Generation Method (CGM) described in Al-Yakoob and Sherali (2013) for solving $\overline{\mathrm{VM}}$, where the transshipment facility was allowed therein to be situated at a point selected from a discretized ad hoc list of three potential locations.
-
• CGM1: The Column Generation Method (CGM) described in Al-Yakoob and Sherali (2013) for solving $\overline{\mathrm{VM}}$, where in the present case, the transshipment facility is allowed to be situated at a point selected from the three best-cost potential facility locations from the set BL obtained via Model AF based on various line segments as described above in Section 6.
-
• CGH: The Column Generation Heuristic (CGH) described in Al-Yakoob and Sherali (2013) for solving VM, where the transshipment facility was allowed therein to be situated at a point selected from a discretized ad hoc list of three potential locations.
-
• CGH1: The Column Generation Heuristic (CGH) described in Al-Yakoob and Sherali (2013) for solving VM, where in the present case, the transshipment facility was allowed therein to be situated at a point selected from the three best-cost potential facility locations form the set BL obtained from solutions of Model AF based on various line segments as described above in Section 6.
-
• ${v_{\mathrm{CGM}}}$: The objective function value of Model $\overline{\mathrm{VM}}$ obtained using Procedure CGM.
-
• ${v_{\mathrm{CGH}}}$: The objective function value of Model VM obtained using Procedure CGH.
-
• ${v_{\mathrm{CGM}\mathrm{1}}}$: The objective function value of Model $\overline{\mathrm{VM}}$ obtained using Procedure CGM1.
-
• ${v_{\mathrm{CGH}\mathrm{1}}}$: The objective function value of Model VM obtained using Procedure CGH1, noting that ${v_{\mathrm{CGH}\mathrm{1}}}$ essentially represents the total cost obtained when implementing Algorithm HA.
-
• $\text{opt-gap}({v_{\mathrm{CGH}}})\equiv 100(\frac{{v_{\mathrm{CGH}}}-{v_{\mathrm{CGM}}}}{{v_{\mathrm{CGH}}}})$, which represents an upper bound on the percentage deviation of ${v_{\mathrm{CGH}}}$ from the (unknown) optimal value for VM.
-
• $\text{opt-gap}({v_{\mathrm{CGH}\mathrm{1}}})\equiv 100(\frac{{v_{\mathrm{CGH}\mathrm{1}}}-{v_{\mathrm{CGM}\mathrm{1}}}}{{v_{\mathrm{CGH}\mathrm{1}}}})$, which represents an upper bound on the percentage deviation of ${v_{\mathrm{CGH}\mathrm{1}}}$ from the (unknown) optimal value for VM.
-
• $\% \hspace{2.5pt}\mathrm{Impr}(\mathrm{HA})\equiv 100(\frac{{v_{\mathrm{CGH}}}-{v_{\mathrm{CGH}\mathrm{1}}}}{{v_{\mathrm{CGH}}}})$, which represent the percentage cost reduction obtained using Algorithm HA over a stand-alone application of Heuristic CGH.
-
• RT: Run-time in seconds (sec), where all runs have been made on a Core${^{\mathrm{TM}}}$ i7 Processor, CPU 4.00 GHz computer having 4 GB of RAM, using the commercial package CPLEX (Version 12.0) as the optimization solver, and with coding in Java.
Note that our implementation of CGH in Al-Yakoob and Sherali (2013), we considered three potential judicious locations for each transshipment facility, although this procedure can generally handle an arbitrary number of potential facility locations. Hence, for the sake of comparison and coordination with our implementation in Al-Yakoob and Sherali (2013), we explore in the first stage of HA a number of line segments (as per the specified data) for Model AF in order to determine the set BL of best-cost potential facility locations as defined in the foregoing section, from which we select the best three cost-effective locations as prescribed by Model AF. We then use these three discrete locations within CGH (referred to in this case as CGH1 in this paper) of Al-Yakoob and Sherali (2013) in the second stage of HA for determining the facility location along with cost-effective shipping schedules. Note that this restricted application of HA was selected for the sake of convenience in implementation of, and comparison with, the procedure of Al-Yakoob and Sherali (2013), but naturally, one would expect further improved results than those reported below by implementing the column generation framework of Al-Yakoob and Sherali (2013) to consider all the potential locations within BL.
Table 1
Computational results for solving Models $\overline{\mathrm{VM}}$ and VM respectively using CGM and CGH.
| Test problem | ${v_{\mathrm{CGM}}}$ (million $) | RT for CGM (seconds) | ${v_{\mathrm{CGH}}}$ (million $) | RT for CGH (seconds) | No. of iterations for CGH | Optimality gap % |
| ${\mathrm{P}_{1}}$ | 28.838014 | 0.29 | 28.838026 | 0.32 | 7 | 4.16117E−05 |
| ${\mathrm{P}_{2}}$ | 35.560039 | 0.34 | 35.854048 | 4.17 | 3 | 0.820016195 |
| ${\mathrm{P}_{3}}$ | 42.310069 | 0.58 | 42.977800 | 6.94 | 18 | 1.553664915 |
| ${\mathrm{P}_{4}}$ | 49.060100 | 0.75 | 50.043900 | 7.30 | 10 | 1.965873963 |
| ${\mathrm{P}_{5}}$ | 55.832607 | 1.13 | 56.885800 | 9.30 | 14 | 1.851416346 |
| ${\mathrm{P}_{6}}$ | 67.195080 | 1.82 | 68.781700 | 6.44 | 15 | 2.306747289 |
| ${\mathrm{P}_{7}}$ | 74.012558 | 2.18 | 75.190300 | 9.16 | 16 | 1.566348319 |
| ${\mathrm{P}_{8}}$ | 96.737634 | 5.06 | 98.179100 | 10.77 | 5 | 1.468200462 |
| ${\mathrm{P}_{9}}$ | 110.372773 | 7.07 | 111.886000 | 11.44 | 15 | 1.352472159 |
| ${\mathrm{P}_{10}}$ | 124.007715 | 12.03 | 125.495000 | 16.30 | 10 | 1.185134866 |
| Average | 68.3926589 | 3.125 | 69.4131674 | 8.214 | 11.3 | 1.406991613 |
Table 1 presents results for the column generation heuristic (CGH) of Al-Yakoob and Sherali (2013) using three potential facility locations selected in an ad hoc fashion. The average optimality gap percentage obtained for CGH1 is about 1.406%, with an associated average overall operational cost of $69.413M. The corresponding average CPU time consumed was about 8.214 seconds. These results confirm the robustness and efficiency of CGH. Improved quality solutions were obtained using Algorithm HA as seen from Table 2, where the average overall operational cost achieved by HA (CGH1) reduced that obtained by CGH by 2.028%, which amounts to an average saving of $2.21M over an average time horizon of 120 days.
Table 2
Computational results for solving Models $\overline{\mathrm{VM}}$ and VM respectively using CGM1 and CGH1.
| Test problem | ${v_{\mathrm{CGM}\mathrm{1}}}$ (million $) | RT for CGM1 (seconds) | ${v_{\mathrm{CGH}\mathrm{1}}}$ (million $) | RT for CGH1 (seconds) | No. of iterations for CGH1 | Optimality gap % | % Impr (HA) |
| ${\mathrm{P}_{1}}$ | 27.245018 | 0.25 | 27.245018 | 0.26 | 4 | 0 | 5.523984201 |
| ${\mathrm{P}_{2}}$ | 33.910047 | 0.37 | 33.930020 | 1.42 | 2 | 0.058865 | 5.366278307 |
| ${\mathrm{P}_{3}}$ | 42.310033 | 0.53 | 42.788000 | 2.17 | 7 | 1.117059 | 0.44162335 |
| ${\mathrm{P}_{4}}$ | 49.060085 | 0.73 | 49.784100 | 4.31 | 17 | 1.454310 | 0.519144191 |
| ${\mathrm{P}_{5}}$ | 55.832300 | 1.00 | 55.832577 | 2.17 | 15 | 0.000496 | 1.851469084 |
| ${\mathrm{P}_{6}}$ | 67.195040 | 1.64 | 67.195100 | 5.23 | 9 | 8.817650 | 2.306718211 |
| ${\mathrm{P}_{7}}$ | 74.012600 | 2.34 | 74.012608 | 6.78 | 17 | 9.619982 | 1.566281821 |
| ${\mathrm{P}_{8}}$ | 96.737706 | 5.14 | 97.957700 | 7.23 | 5 | 1.245429 | 0.225506243 |
| ${\mathrm{P}_{9}}$ | 110.372593 | 7.46 | 110.373000 | 9.63 | 9 | 0.000368 | 1.352269274 |
| ${\mathrm{P}_{10}}$ | 124.077994 | 11.54 | 124.080000 | 18.11 | 13 | 0.000162 | 1.127534962 |
| Average | 68.0753416 | 3.1 | 68.3198123 | 5.731 | 9.8 | 2.2314321 | 2.028080964 |
As a point of interest, we mention that the average optimality gap percentage obtained when solving AF using RHA by itself (based on the best location within BL) turned out to be 15.21%, with an associated average overall operational cost of $78.59M. The corresponding average CPU time required was 758.09 seconds. The average cost reduction in this first stage solution value obtained by applying the second stage of the proposed hybrid heuristic HA was given by 13.08%, thus underlining the importance of applying CGH of Al-Yakoob and Sherali (2013) in order to select among a judicious set of potential facility locations and to determine associated cost-effective shipping schedules. Naturally, by implementing CGH1 to select among all the potential locations in the set BL as prescribed by Model AF within the context of the hybrid heuristic HA, as opposed to just three locations as used in our implementation, would likely further benefit the quality of the derived solution.
8 Summary, Conclusions, and Future Research
This paper explores mathematical modelling approaches for determining an optimal location for a transshipment facility in a single source-destination vessel scheduling and transportation-inventory problem. The problem is concerned with transporting a product from a source to a destination based on a stream of consumption rates at a client delivery destination. Different cost components are taken into consideration pertaining to the daily operation of vessels, chartering expenses, and penalties associated with undesirable storage levels at destination. This research effort is a continuation of the authors’ work in Sherali and Al-Yakoob (2006a, 2006b) and Al-Yakoob and Sherali (2013), which focuses on exploring potential savings that can be achieved by strategically locating a transshipment facility for use during operation. In particular, this research extends the foregoing work by incorporating modelling considerations related to determining a location for a transshipment facility within a restricted continuous region for the single source-destination operation. A two-stage hybrid approach (HA) has been proposed in this paper, which first determines a set of cost-effective facility locations and then uses these locations in CGH1 to further achieve cost reduction. Based on our 10 test problems, improved overall operational costs were achieved using HA over a stand-alone application of CGH as in Al-Yakoob and Sherali (2013), where the average overall operational cost obtained by the former method was reduced by 2.028% in comparison with that obtained via the latter, amounting to an average saving of $2.21M over an average time horizon of 120 days.
As an extension to this research effort, we recommend exploring for both single as well as multiple source-destination problems alternative models that can be directly solved using a column generation approach in lieu of using the proposed two-stage hybrid algorithm.
Authors’ contribution. S. Al-Yakoob: Manuscript Writing, Literature Search and Review, Industrial Consultation related to the studied problem, Mathematical Modelling and Analysis, and Computational Results. H.D. Sherali: Manuscript Writing and Mathematical Modelling and Analysis.
