Abstract

This paper presents a novel application of operations research to support decision making in blood distribution management. The rapid and dynamic increasing demand, criticality of the product, storage, handling, and distribution requirements, and the different geographical locations of hospitals and medical centers have made blood distribution a complex and important problem. In this study, a real blood distribution problem containing 24 hospitals was tackled by the authors, and an exact approach was presented. The objective of the problem is to distribute blood and its products among hospitals and medical centers such that the total waiting time of those requiring the product is minimized. Following the exact solution, a hybrid heuristic algorithm is proposed. Computational experiments showed the optimal solutions could be obtained for medium size instances, while for larger instances the proposed hybrid heuristic is very competitive.

1. Introduction

This paper presents a novel application of operations research to support decision making in blood distribution management with the objective of minimizing the total waiting time of those hospitals and medical centers requesting for blood products. The increasing demand for healthcare services coupled with their importance and higher costs make it obligatory to better utilize the medical resources and facilities. The rapid and dynamic increasing demand for blood, criticality of the product, storage, handling, and distribution requirements and limitations, and the different geographical locations of hospitals and medical centers have made blood distribution a complex and very important problem. Furthermore, increasing rate of surgeries coupled with the new advances in healthcare have magnified the complexity and importance of an efficient blood distribution system. If the right blood products are not available at the hospital and medical centers at the right time, then health issues or delays of operations may arise, which results in extra days of hospitalization and cost. Besides, overstocking blood at hospitals and medical centers leads to low utilization, as most blood products may only be used to a patient of the same blood type within 21 days of collection. For these reasons, blood must be collected regularly. Thus, any low utilization increases costs and wastes the scarce blood resource. The latter may have fatal consequences.

The human blood is one of the most important components of a healthcare system. According to [1], the most important reasons for the need for this vital medical resource are as follows.(i)There is no exact substitute for human blood.(ii)Blood and products made from blood play an important role in advances in lifesaving techniques.(iii)Blood products cannot be stored for an indefinite period (e.g., according to [2], red blood cells must be used within 35–42 days of collection; platelets have shorter life, as they must be used within five days of collection).

Blood is composed of many components (red cells, white cells, platelets, plasma). Each of these components supplies a separate function in the human organism and has a different use in medical treatment. Despite the various requirements to the number of unit of blood, almost all critical medical treatments and operations require this vital life resource, among which are accidental victims and injuries, surgeries, organ transplantations, several cancer treatments, and so forth. For instance, 8 units of platelets may be required on a daily basis by a patient undergoing leukemia cancer treatment [1].

The previous studies on the blood products distribution have focused on the inventory related problems [37]. However, there is a lack of study regarding the problem of distributing blood and blood products among hospitals and medical centers with the objective of minimizing the total waiting time of hospitals waiting for blood arrival. Delen et al. studied blood supply chain management by analyzing inventory consumption patterns and supply chain status. They implemented the approach in the supply chain facilities at the two United States Air Force Bases [1]. Katsaliaki performed a simulation study towards a more cost-effective management of blood supply chain in the United Kingdom. The study revealed substantial improvements in inventory and distribution operations involved and resulted in cost reductions and increased safety [8]. Katsaliaki and Brailsford studied ordering policies with respect to reductions in shortages and wastage, associated costs, as well as improving service levels and safety by applying a discrete-event simulation model [9]. Hemmelmayr et al. studied the problem of a cost-effective delivery of blood products to Austrian hospitals. They presented solution approaches based on integer programming and Variable Neighborhood Search meta-heuristic [10]. Although less related, Rego and Sousa performed a thorough study on the design of hospital supply chains systems and developed a hybrid Tabu Search-Variable Neighborhood Search meta-heuristic for the problem [11].

In this study, a real-blood distribution problem was considered. The case was experienced by authors in the city of Tehran, Iran. The authors were supposed to derive practical solutions, preferably optimal distribution of blood products among hospitals and medical centers, with the objective of minimizing total waiting time of those requiring the product. As solution procedures, a mixed-integer programming formulation was presented followed by a hybrid heuristic to find optimal distribution of blood products among hospitals and medical centers such that the total waiting time of those requiring the product is minimized. Computational experiments showed the mixed-integer programming formulation is quite capable of finding optimal solutions for medium-size instances, and for the real cases where 24 hospitals were serviced, while for larger random instances, the hybrid heuristic derives near optimal solutions very quickly. The remainder of this paper is organized as follows. In Section 2, we give an exact definition of the problem together with the notations used in this paper. In Section 3, we develop a new mixed-integer programming formulation for the problem. Following incapability of the exact procedures in deriving optimal solution for large instances in reasonable time, Section 4 is devoted to the hybrid heuristic algorithm developed to solve the large-scale instances of the problem. In Section 5, computational experiments are reported. The paper ends with conclusion.

2. Problem Definition

This section defines the problem, and the notations used throughout this paper. Apart from many improvements applied to blood donor service processes since its foundation in 1921, blood products supply chain and distribution process are composed of the four steps:(1)collection of blood from donors,(2)testing, processing and storage of blood,(3)delivery of blood products,(4)consumption of blood products to recipients.

The problem studied in this paper is inherent in process three, and it is to deliver the blood to hospitals and medical centers upon their requests such that the total waiting time of those requesting for is minimized. For the sake of simplicity, when “hospital” is used, we mean both hospitals and medical centers, or any other medical institute where blood and blood products are required. Typically, a central blood transfusion depot distributes blood among hospitals. Apart from impartiality in blood distribution to hospitals when a central depot is implemented, the total costs and waiting times can be minimized while managing the whole operation is much simpler. Besides, it is not appropriate to keep the blood for a long time and several blood products have a very short life, hence, keeping them for even days results in unusable blood products which cannot be transferred to patients. These limit banking large volume of blood by hospitals. Apart from these limitations which build this problem very interesting and challengeable, the operational limitations avoid a single blood distribution vehicle to cover all hospitals on a single route. At the most, these operational limitations arise from the fact that the blood should be delivered to hospitals before the latest time (an upper bound on the delivery time), and also each vehicle has a limited capacity. In fact, in reality each vehicle covers a set of hospitals where a route is built.

An example of the problem is illustrated in Figure 1 where a set of 50 hospitals are served by two vehicles. The vehicles leave from blood transfusion depot (the black rectangle). Initial studies show that at least two vehicles are required to cover these 50 hospitals. Thus, two decisions should be made: Which hospitals should be assigned to each vehicle? What the objective function would be?, and How the permutation of visiting each hospital by each vehicle (in fact, in each route) should be determined? (Again, what the objective function would be?). In this study, these two decisions are made separately by employing two different mathematical formulations, where, at first, hospitals are assigned to blood distribution vehicles (or simply constructing routes), and then permutation of visiting hospitals by each vehicle are optimized, according to the appropriate objective function criterion. To better understand the second problem, that is, optimal permutation of visiting hospitals in each route, the following example is provided.

Assume the problem of Figure 1(b), where a set of 15 hospitals is covered by a single vehicle (vehicle 1). Without loss of generality and for the sake of better understanding, assume that the travel time is our only component of the objective function. Table 1 shows these travel times as a symmetric travel time matrix. We hold “0” as the blood transfusion depot of vehicles, or simply depot according to the traveling salesman problem (TSP) and vehicle routing problem (VRP) literatures. We define the “arrival time” to a hospital as the time required for the vehicle to visit this hospital from depot for the first time. For instance, according to Figure 2(a), the permutation of visiting hospitals is 𝑝3𝑎={0,44,16,36,39,37,25,47,41,23,27,35,45,48,7,17,0}.(2.1) Thus arrival time to hospital 44 would be 64 (see Table 1). This is the total waiting time of hospital 44 to be visited by the vehicle 1.

Similarly, the arrival time to the hospital 16 would be the arrival time to hospital 44 (as hospital 16 cannot be visited any sooner than hospital 44 in permutation 𝑝3𝑎) plus the travel time from hospital 44 to hospital 16. Thus, 𝑎16=44+𝑡44,16=64+7=71.

Obviously, this arrival time concept is similar to the waiting time concept. In fact, what we are looking for is to decide how to locate the hospitals so as to minimize the total waiting time of all hospitals assigned to a vehicle. It is not surprising that the most important issue when developing healthcare systems, especially in an emergency situation, is the waiting time. Figures 2(a) and 2(b) illustrate two examples of waiting time issue, where changing only direction of visits in Figure 2(b) in comparison with Figure 2(a) leads to another solution with a different objective function value (assuming routes are constructed). Minimizing the total waiting time of hospitals waiting for blood can be modeled within the framework of the traveling repairman problem ((TRP), (other names are minimum latency problem (MLP) and travelling delivery-man problem (TDP)). The TRP has been proved to be NP-Hard and even more difficult than the TSP [12, 13]. Recently, Salehipour and Sepehri have developed a strong mixed-integer mathematical programming formulation for the problem based on the flow network problem [14]. Fischetti et al. [15] and Eijl [16] have developed two mathematical programming formulations for the TRP, however, their formulations are much weaker than that of [14], as their models cannot solve even small instance size of 15 hospitals in a reasonable time. Sarubbi and Luna presented a different mixed-integer mathematical programming model at the International Network Optimization Conference in 2007 [17]. However, their model still lacked certain sets of constraints and resulted in infeasibility. Mendez-Diaz et al. [18] formulated TRP based on the linear ordering problem. Although model’s performance is much better than previous ones, its computational time is not comparable to that of Salehipour and Sepehri [14]. Later, we report on dimensions (total number of decision variables and constraints), and the computational performance of these models.

Note that here a solution is represented as a permutation of hospitals starting from 0 and ending in 0, that is, 𝑝={0,,𝑖,𝑗,,0}. From this representation, the objective function value can be easily calculated. Assume the solution 𝑝={0,1,2,,𝑛,0} with 1 being the first hospital in the solution; 2, the second, and so on. The total arrival time (equivalently waiting time) for this solution can be calculated as 𝐹𝑝=𝑛𝑖=1𝑎𝑖,(2.2) where 𝑎𝑖=𝑗𝑖𝑗=1𝑡𝑗1,𝑗.(2.3) Combining (2.2) and (2.3) and simplifying it result 𝐹𝑝=𝑛𝑖=1((𝑛+1)𝑖+1)𝑡𝑖1,𝑖+𝑡𝑛,0.(2.4) Equation (2.4) offers the advantage of making the calculations much easier and quicker, especially in heuristic algorithms. Let us back to Example of Figure 2. The permutation 𝑝3𝑎 has a better (lower) objective function value than the permutation 𝑝3𝑏: 𝐹𝑝3𝑎=𝑛𝑖=0𝑎𝑖𝐹=16(64)+15(7)+14(14)+13(8)+12(5)+11(24)+10(6)+9(7)+8(14)+7(44)+6(25)+5(13)+4(35)+3(11)+2(35)+1(66)=2820,𝑝3𝑏=𝑛𝑖=0𝑎𝑖=16(66)+15(35)+14(11)+13(35)+12(13)+11(25)+10(44)+9(14)+8(7)+7(6)+6(24)+5(5)+4(8)+3(14)+2(7)+1(64)=3606.(2.5)

2.1. Notations

We assume that the time component consists of the total traveling time of blood distribution vehicles between every two arbitrary hospitals, 𝑖 and 𝑗, and that the service time (the delivery time) is assumed to be same for all vehicles and hospitals, and hence is ignored. Even not ignored, this delivery time can be dealt with easily [19]. We also assume “0” as the initial position of all vehicles.

2.1.1. Decision Variables

The following hold.𝑦𝑖𝑗 is a decision variable which takes 1 if the vehicle travels from hospital 𝑖 to hospital 𝑗; and 0, otherwise.𝑥𝑖𝑗 is a decision variable which takes 𝑛𝑢+1, if hospital 𝑗 is visited after hospital 𝑖 and is visited by a vehicle in order 𝑢; and 0, otherwise.𝑝𝑢𝑗 is a decision variable which takes 1 if hospital 𝑗 is visited in order 𝑢; and 0, otherwise.𝛿𝑣𝑖 is a decision variable which takes 1 if hospital 𝑖 is visited by vehicle 𝑣; and 0, otherwise.𝛾𝑣 is a decision variable which takes 1 if vehicle 𝑣 is used; and 0, otherwise.

2.1.2. Parameters and Indices

The following hold.𝑖 Hospital 𝑖, 𝑖=0,,𝑛, where 𝑖=0 is a dummy hospital and corresponds to the blood transfusion depot of vehicles (vehicles’ home).𝑛+1 is total number of hospitals (including the depot).𝑎𝑖 is arrival time to hospital 𝑖, where 𝑎0=0.𝐹𝑝 is total cumulative travel time for route 𝑝 (permutation 𝑝).𝑐𝑣 is the cost of using vehicle 𝑣, 𝑣=1,,𝑚.𝑈𝑇 is the upper bound on the total travel time of each vehicle (since we have 𝑚 vehicles, thus we have 𝑚 routes).𝑡𝑖𝑗 is a parameter which states the total travel time from hospital 𝑖 to hospital 𝑗.

3. Problem Formulation

This section describes one of the two major contributions of this study by proposing an exact solution approach for the problem. The idea is to decompose the problem into subproblems using set covering formulations, and then to solve optimally each subproblem using TRP formulations. Thus, at first hospitals are assigned to vehicles (route construction), and then visiting order of hospitals inside each route is optimized. Computational complexity of the problem restricts developing a single combined mathematical programming formulation which assigns hospitals to vehicles while minimizes total waiting time inside each route.

3.1. Constructing Routes: A Location Covering Model

The assignment of hospitals to vehicles can be accomplished using the set covering formulations. One of such formulations has appeared in [20] where the goal is to locate a least-cost set of facilities such that each customer (hospital) can be reached within a maximum allowed travel time from the closest facility (blood transfusion depot). To apply this formulation to construct blood distribution routes, however, small modifications are required. These modifications are performed by imposing a set of constraints inspired from operational limitations and problem’s nature (see constraints (3.4) below), and several changes into the model’s parameters. Furthermore, as all vehicles leave from a central blood transfusion depot, the concept of closest facility as appeared in the location covering formulation is not applicable here. However, applying the constraints (3.4) coupled with the location covering model provided in [20] (see Model 1 below), ensures assignment of a set of hospitals to each vehicle where minimum travel time is satisfied.

Model 1
Minimize𝑧=𝑚𝑣=1𝑐𝑣𝛾𝑣+𝑛𝑚𝑖=0𝑣=1𝑡0𝑖𝛿𝑣𝑖(3.1) subject to 𝑚𝑣=1𝛿𝑣𝑖1,𝑖=0,,𝑛,(3.2)𝑛𝑖=0𝛿𝑣𝑖||||𝛾𝑛+1𝑣,𝑣=1,,𝑚,(3.3)𝑛𝑖=0𝑡0𝑖𝛿𝑣𝑖𝛿𝑈𝑇,𝑣=1,,𝑚,(3.4)𝑣𝑖,𝛾𝑣{0,1},𝑖=0,,𝑛,𝑣=1,,𝑚.(3.5) Constraints (3.2) guarantee that all hospitals are visited, while constraints (3.3) ensure that if vehicle 𝑣 is not set up, then no hospital can be covered by it. Constraints (3.4) ensure that total travel time of each route associated with each vehicle satisfies operational limitations. Constraints (3.5) state that all 𝛿𝑣𝑖 and 𝛾𝑣 are binary.

3.2. Visiting Hospitals: A Traveling Repairman Model

After assigning hospitals to vehicles and constructing routes, optimal permutation of visiting hospitals with respect to minimizing total waiting time (equivalently total cumulative traveling time) of hospitals should be determined. This section describes a novel mixed-integer programming formulation for this purpose (Model Flow-SP). The proposed formulation is strong enough and allows optimal solutions for instances of up to 25 hospitals in a single route, and very close to optimality for instances with 30 hospitals.

Model Flow-SP
Minimize𝑧=𝑛𝑛𝑖=0𝑗=0𝑡𝑖𝑗𝑥𝑖𝑗(3.6) subject to 𝑛𝑖=0𝑦𝑖𝑗=1,𝑗=1,,𝑛,(3.7)𝑛𝑗=0𝑦𝑖𝑗=1,𝑖=1,,𝑛,(3.8)𝑛𝑖=1𝑦𝑖0=1,(3.9)𝑛𝑖=1𝑦0𝑗=1,(3.10)𝑛𝑗=0𝑥0𝑗=𝑛+1,(3.11)𝑛𝑗=0𝑥𝑗0=1,(3.12)𝑛𝑖=0𝑥𝑖0𝑛𝑗=0𝑥0𝑗=1(𝑛+1),(3.13)𝑛𝑖=0𝑥𝑖𝑘𝑛𝑗=0𝑥𝑘𝑗=1,𝑘=1,,𝑛,(3.14)𝑛𝑛𝑖=0𝑗=0𝑥𝑖𝑗=𝑛+1𝑠=1𝑥𝑠,(3.15)𝑖𝑗(𝑛+1)𝑦𝑖𝑗𝑦,𝑖,𝑗=0,,𝑛,(3.16)𝑖𝑗𝑥{0,1},𝑖,𝑗=0,,𝑛,(3.17)𝑖𝑗0,𝑖,𝑗=0,,𝑛,(3.18)𝑛𝑢=0𝑝𝑢𝑗=1,𝑗=0,,𝑛,(3.19)𝑛𝑗=0𝑝𝑢𝑗𝑝=1,𝑢=0,,𝑛,(3.20)00=1,(3.21)𝑛𝑗=0𝑝𝑛𝑗=1,(3.22)𝑛𝑘=0𝑥𝑗𝑘(𝑛+1(𝑖+1))𝑝𝑖𝑗0,𝑖,𝑗=0,,𝑛,(3.23)𝑛𝑗=0𝑦𝑖𝑗𝑛𝑢=0𝑝𝑢𝑖𝑝0,𝑖=0,,𝑛,(3.24)𝑢𝑖{0,1},𝑢,𝑖=0,,𝑛.(3.25) Constraints (3.7) and (3.8) are the assignment problem constraints that ensure every hospital is visited exactly once. Constraints (3.9) and (3.10) ensure that the vehicle starts from depot and finishes at depot. Redundant constraints (3.11) and (3.12) together with constraints (3.13), (3.14), and (3.15) define a network flow problem to remove subtours, that is, ensuring feasibility of solutions. Constraints (3.19) to (3.25) are scheduling constraints working according to the visiting order of the hospitals. Constraints (3.19) and (3.20) ensure that there is just one hospital visited in position 𝑖 and vice verse. Constraints (3.21) and (3.22) ensure that the vehicle starts from depot and goes back to depot at the end of its service. Constraint (3.16), and constraints (3.23) and (3.24), links variables 𝑥𝑖𝑗 to variables 𝑦𝑖𝑗, and variables 𝑦𝑖𝑗 to variables 𝑝𝑢𝑖, respectively. The fact that all 𝑥𝑖𝑗, 𝑝𝑢𝑖, and 𝑦𝑖𝑗 variables are binary and nonnegative, respectively, is assured by constraints (3.17), (3.25), and (3.18).
Table 2 reports the dimensions of four mathematical models developed for TRP, namely Model of Fischetti et al. [15], Model of Eijl [16], Model of Mendez-Diaz et al. [18], and that of Salehipour and Sepehri [14] (also Model Flow-SP).

Apart from Table 2, another comparison regarding each model’s performance is reported in Table 3. In this table, 50 problems in 10 different sizes ranging from 10 nodes to 50 nodes were solved. In each size, 5 instances were considered, where nodes’ coordinates were randomly generated, and Euclidean costs were calculated and rounded down to the nearest integer. For each size, mean, minimum, and maximum computational time in seconds (over the five instances), and mean, minimum, and maximum gap in percent were reported. Note that the limit on the computational time is set to be 1 hour. The commercial solver Cplex 9.0 from ILOG was implemented.

According to the table, the best results were reported by Model Flow-SP developed by Salehipour and Sepehri [14]. The performance of Model of Mendez-Diaz et al. is quite promising. Even ignoring the computational time, and gaps, Model of Mendez-Diaz et al. is not still comparable to Model Flow-SP, as the largest size solved by Model Flow-SP has 25 nodes.

4. Heuristic Model

Despite a strong mixed-integer formulation developed in Section 3.2, the problem of finding optimal permutation hospitals in a route (traveling repairman problem) is still very difficult to solve. To derive near-optimal solutions for medium and larger instances, an efficient hybrid heuristic algorithm is developed in this section. The basic idea is to control the improvement procedure accomplished by the variable neighborhood search (VNS) algorithm using the simulated annealing (SA) algorithm. The basis of the heuristics developed in the following sections require finding the closest, second-closest, and so forth hospital to any given hospital. We, therefore, preprocess the travel time data and maintain in memory a proximity matrix that contains one row per hospital. The 𝑖th row of the matrix contains all hospitals (except for hospital 𝑖) in order of proximity to hospital 𝑖. Followings are through details of this hybrid algorithm.

4.1. The Construction Algorithm

The construction heuristic is the greedy randomized adaptive search procedure (GRASP) developed by Feo and Resende [21, 22]. The motivation for using the GRASP for the TRP is the observation that the first few hospitals in a route are more important than the later ones. This immediately follows from (2.4), from which it is clear that the 𝑖th distance in the solution is multiplied by a factor (𝑛+1)𝑖+1. Hence, it is expected that a greedy algorithm performs well. A completely greedy algorithm on the other hand can be expected to miss some interesting opportunities.

As the most distinguishing characteristic of the GRASP where greediness is combined with randomness in the construction phase, a restricted candidate list (RCL) is built by selecting a subset of all elements (hospitals) in a greedy fashion. Assuming a minimization problem, the RCL contains the elements whose incorporation into the partially built solution would yield the smallest increase in the objective function value. From the RCL, an element is then selected at random, after which the RCL is updated to reflect the fact that a new element was added to the solution and is no longer available for selection. Selection of an element and update of the RCL are repeated until a complete solution has been built. From this solution, an improvement procedure starts until a local optimum is found. The size of the RCL, 𝛼, is a parameter of the GRASP algorithm that controls the balance between greediness and randomness. If 𝛼 is small, the search is relatively greedy. If 𝛼 is large, it is relatively random.

The RCL concept can be easily translated to the TRP, and efficiently implemented using the proximity matrix. Starting from 0, the RCL is filled with the 𝛼 closest hospitals. From this RCL, a random hospital is chosen, and then the RCL is filled with the 𝛼 hospitals closest to it, removing any hospitals that have already been selected. If less than 𝛼 hospitals remain, then the size of the RCL is decreased. Algorithm 1 outlines the GRASP construction algorithm for the TRP.

𝑈 { 0 , 𝑣 1 , 𝑣 2 , , 𝑣 𝑛 } ;
𝑣 𝑐 0 ;
Repeat
Create RCL with α vertices 𝑣 𝑖 𝑈 closest to 𝑣 𝑐 ;
Select random vertex 𝑣 𝑟 RCL;
𝑈 𝑈 { 𝑣 𝑟 } ;
𝑣 𝑐 𝑣 𝑟
Until   𝑈 = ;

4.2. The Improvement Procedure: SA + VNS Algorithms

The improvement procedure includes the five local searches, acquired from the routing literature and modified for this problem, changed in a systematic way using the variable neighborhood search (VNS) meta-heuristic [2326]. This VNS algorithm is itself controlled by the simulated annealing (SA) meta-heuristic to avoid being trapped in local optima.

The VNS meta-heuristic systematically explores different neighborhood structures. The main idea underlying the VNS is that a local optimum relative to a certain neighborhood structure is not necessarily a local optimum relative to another neighborhood structure. For this reason, escaping from a local optimum can be accomplished by changing the neighborhood structure. Although this is not required, many implementations of VNS use a sequence of nested neighborhoods, 𝑁1 to 𝑁𝑘max, in which each neighborhood in the sequence is a superset of its predecessor, that is, 𝑁𝑘𝑁𝑘+1. The neighborhoods used in our VNS do not possess this property. Furthermore, usually in the VNS algorithm larger neighborhoods are only examined when all smaller neighborhoods have been depleted, that is, when the current solution is a local optimum of all smaller neighborhoods. Another option is to use the neighborhoods in the order of their effectiveness with respect to the problem at hand, which can be established by a small pilot study. Algorithm 2 presents the pseudocode for a basic VNS.

Input: Initial solution 𝑥 ;
Until the stopping criterion is met do
 Initialize: 𝑘 1 ;
Until   𝑘 = 𝑘 m a x   repeat
  Shake: Generate ́ 𝑥 randomly from 𝑁 𝑘 ( 𝑥 ) and set it as the incumbent solution;
  Local search: Apply some local search method on ́ 𝑥 ; Denote with 𝑥 the so obtained local
  optimum;
  If   𝑥 is better than the 𝑥 then
   𝑥 𝑥 and set 𝑘 1 ;
  Else
   𝑘 𝑘 + 1 (switch to another neighborhood);

Based on our previous experience, we set 𝑘max=5 in the hybrid algorithm. These five neighborhoods were applied according to the order presented in Table 4.

Introduced by Kirkpatrick et al. [27], the SA algorithm is one of the well-known meta-heuristic algorithms where trapping in local optima is avoided by sometimes accepting a neighborhood move which worsens the value of the objective function. The acceptance or rejection of a move (even the worse moves) is determined by a sequence of random numbers, but with a controlled probability. The probability of accepting a move, causing an increase Δ in the objective function 𝑓, is called the acceptance function and is normally set to 𝑒(Δ)/𝑇, where 𝑇 is a control parameter corresponding to the temperature in the analogy with physical annealing. This acceptance function implies that small increases in 𝑓 are more likely to be accepted than large increases. When 𝑇 is high, most moves will be accepted, but as 𝑇 approaches zero most uphill moves will be rejected. Usually, the SA starts with a relatively high value of 𝑇 to avoid being trapped in a local optimum too early. The algorithm proceeds by attempting a certain number of neighborhood moves at each temperature (here, in a VNS algorithm), while the temperature gradually drops. We implemented the SA algorithm as a controlling scheme on local searches to avoid being trapped in local optima. The outline of the SA is brought in Algorithm 3, and the outline of the hybrid algorithm is shown in Algorithm 4.

Input: Initial solution 𝑥
Select an initial temperature 𝑇 > 0 ;
Repeat
Local search: ́ 𝑥 a r g m i n 𝑁 ( 𝑥 ) ;
Calculate Δ = 𝑓 ( ́ 𝑥 ) 𝑓 ( 𝑥 ) ;
If Δ < 0 then 𝑥 ́ 𝑥 ;
Elseif random(0,1) < 𝑒 Δ / 𝑇 , then 𝑥 ́ 𝑥 ;
Else
Reject ́ 𝑥 ;
𝑇 = 𝑇 ( 𝑐 ) ;
Until stopping condition is met;

Initial procedure GRASP, report the best found solution 𝑥 ;
Initial procedure SA + VNS;
Repeat
Local search: apply VNS procedure with the five neighborhoods, report the incumbent
solution ́ 𝑥 ;
Calculate Δ ;
If Δ < 0 then 𝑥 ́ 𝑥 ;
Elseif random(0,1) < 𝑒 Δ / 𝑇 , then 𝑥 ́ 𝑥 ;
Else Reject ́ 𝑥 , and modify temperature;
Until stopping criterion is met;

4.3. Neighborhoods

Our search uses five neighborhood structures: swap-adjacent, swap, remove-insert, 2-opt, and 3-opt. The swap heuristic attempts to swap the positions of each pair of hospitals in the permutation. The remove-insert heuristic examines randomly each hospital 𝑖 in the solution and places it at the end of the permutation. The swap-adjacent heuristic, a subset of the swap heuristic, attempts to swap each pair of adjacent hospitals in the permutation. The 2-opt heuristic removes each pair of edges (partial routes between two hospitals) from the solution and reconnects the nodes (hospitals). The 3-opt heuristic removes each triplet of edges from the solution and reconnects the solution in such a way that the orientation of the solution parts is preserved.

4.4. Neighborhood Reduction

Several attempts are made to decrease the size of the neighborhoods to examine during the improvement phase. The implemented neighborhood reduction scheme is based on the observation that in a reasonably good solution, most improving moves will involve hospitals that appear in relative proximity to each other in this solution. For example, in a good solution it is unlikely that a swap will yield a better solution if it exchanges a hospital that appears close to the depot with one that appears far from it. We, therefore, introduce a proximity factor 𝛽 to determine the maximum distance between hospitals that may be involved in a move. The factor 𝛽 has a slightly different influence, depending on the move. The neighborhood reduction strategy is not used for the swap-adjacent and remove-insert heuristics.

5. Computational Results

In this Section, we provide detailed computational experiments of the mathematical formulation provided in Section 3, and the hybrid heuristic of Section 4. In this study, two datasets were considered, random datasets, and a real dataset. All computational experiments were carried out on a Pentium 4 PC with 2 GB of memory and 2.0 GHz of CPU. We employed solver Cplex 9.0 from ILOG. The hybrid heuristic was coded in C++.

5.1. Random Datasets

A set of 30 random instances were generated with sizes ranging of 50, 100, and 150 hospitals (for each size, 10 problems were generated). For each problem, travel time among hospitals has been generated randomly using uniform distribution in [1,100]. Table 5 shows the general details of the instances generated. In each problem, a maximum number of 10 blood distribution vehicles are considered, although, according to Model 1 not all the vehicles may be implemented in optimal solution. In all instance, we set 𝑈𝑇 (upper bound on the total travel time of each vehicle) to be 1500.

To tune the parameters of our hybrid algorithms, that is, the GRASP algorithm parameter, and the SA algorithm parameters, we refer the interested readers to [14, 28] where a comprehensive study regarding this tuning process for the TRP has been accomplished. The best values for these parameters are shown in Table 6.

Table 7 reports the complete computational results of both mixed-integer mathematical programming and the hybrid heuristic. The forth column of the table, “OFV of Model 1” shows the objective function value associated with Model 1, that is, the minimum cost associated with the constructed routes. The column “OFV of Model Flow-SP” refers to the computational results of the Model Flow-SP, that is, optimal permutation of visiting hospitals in each route. All reported computational times are in seconds, and we set the upper bound on the CPU time to be 500 seconds. For each instance, under the “OFV of Model Flow-SP,” we have shown average, minimum, and maximum gap and time associated with each of the solutions methods (mixed-integer programming and hybrid heuristic). All reported gaps for the exact solution are those from Cplex solver. This gap is calculated from the lower bound derived by the solver over the best objective function found by the solver (obviously, for the optimal solution this gap is 0%). In fact, Cplex solver proves optimality of the solution by using this gap. Thus, this gap can be a measure of strength of the formulation as the convergence rate of this gap towards 0% in a reasonable amount of computational time is of importance. It is worth mentioning that Model Flow-SP yielded an average gap of 14.63% for instances with 50 hospitals, and around 18%, for instances with 100 and 150 hospitals, respectively. The average of “Min. Gap” is even much lower as it is below 8%.

According to the table, the hybrid heuristic exhibited much better performance, both in CPU time and in reported gap, especially when the size of problem increased. This reveals the strength of the developed heuristic for the problem. Here, the reported gap was calculated from the best lower bound found by Cplex solver over the best objective function found by the heuristic. Furthermore, maintaining the same quality, even in many cases with a lower gap, the hybrid heuristic rarely reached the upper bound on the computational time while for instances with 100 and 150 hospitals, the mixed-integer programming formulation reached this upper bound in several cases. Note that according to the table, Model Flow-SP is very strong as for instances with 50 hospitals, its performance is well enough and its average gap is below 15%. But due to the complexity of the problem, still the problems of 50 hospitals are not solvable in a reasonable time.

The developed hybrid-heuristic performed quite well, and tracking back its performance, we observed that implementing the GRASP has been a very good initialization. Also powerful neighborhoods in the VNS meta-heuristic coupled with the SA algorithm has magnified this performance which resulted in near optimal solutions very quickly.

5.2. Real Dataset: A Case Study in Tehran, Iran

Here, the method of Section 3 is applied to a real blood distribution problem, experienced by the authors. In this study, 24 hospitals were considered in the city of Tehran, capital of Iran. Generally, in Tehran a single blood transfusion depot supplies the demands of all hospitals, in a way that each hospital sends its vehicle to deliver the demand. The data of this study includes the travelling distances to reach each hospital, the daily demand (blood products) of each hospital, the maximum travelling distance allowed for vehicles, and the vehicles’ capacity. Note that when considering the travelling cost of each vehicle, we take the travelling distance instead, usually by taking into account the approximate travelling cost per kilometer.

The travelling distances, in kilometer, are given in the appendix. For the sake of simplicity, we assumed the table is symmetric. The distances were extracted by Google Map service of Google. Typically, this service offers several routes between every two points, among which we chose the shortest one.

We were restricted to provide the daily demands of hospitals to public domain. As for vehicles’ capacity, we set 600 units of blood products per vehicle, although different vehicles may hold different capacities. On a standard vehicle, it is common in practice that the maximum travelling distance allowed for vehicle is reached before the capacity. Hence, from practical viewpoint, the capacity may less affect the final decision. It is worth mentioning that the maximum travelling distance allowed differs when different numbers of vehicles are hired. In fact, as the demands of all hospitals should be met, this upper bound on travelling distance should be set such that all demands are met. The maximum travelling distance allowed when hiring different numbers of vehicles are shown in Table 8. Except the case when there is only one vehicle, the data of the table reflect the concern of several hospitals.

According to the number of hospitals considered in this study, it is expected that at most three vehicles could satisfy the demands of all hospitals. The routes, constructed by hiring different numbers of vehicles, are depicted in Figures 3, 4, and 5. Table 9 elaborates more on the results.

In the table, the first column shows the number of vehicles hired for supplying demands while the second column shows the objective function of subproblem of Section 3.1. The remaining columns which are same for each vehicle are interpreted as objective function of the TRP (column OFV), total distances travelled by vehicle in kilometers (column Trv. Dis.), and total travelling time of vehicle ignoring the delivery time (column Trv. Time).

5.3. Current Situation

Currently, in Tehran, each hospital maintains its own vehicle to deliver its demand. However, an investigation performed by the authors revealed that the vehicle is not properly suited for this purpose, which may affect the quality of the service. The disadvantages of this scenario (Scenario 2) are the follows.(i)The central blood transfusion of Tehran which distributes blood products to all hospitals is located in a highly crowded part of the city, where the traffic is limited. Furthermore, only a few hospitals hire vehicles allocated for blood products delivery, which as emergency vehicles, are also allowed to freely enter this area. For hospitals not hiring such vehicles, the most appropriate vehicle would be taxis. Although, taxis are allowed to enter this area, they are not considered as emergency vehicles. Thus, in case of an accident, traffic jam, or any other similar circumstances, they will end up in delay, even to several hours, as an unofficial report by one hospital states.(ii)When hospitals hire taxis, obviously, the vehicles are not specifically designed for the purpose of blood products distribution. Furthermore, the driver has not been trained regarding any special conditions of keeping and transporting blood products (e.g., appropriate temperature).(iii)The human errors could make the situation much worse. For instance, the driver cannot distinguish between different blood products, which probably results in wrong delivery.

Apart from these disadvantages, this scenario incurs higher costs compared to case where specific vehicles are sent to hospitals from blood transfusion depot (Scenario 1, see also Section 2). For this purpose, Table 10 illustrates the associated costs.

The costs of the table are calculated by twice of the travelling distances in kilometers plus 10% of this cost, as for delays in delivery. As the table shows, the total cost is 450.45 (kilometers). Interpreting this cost as the total travelling times of all vehicles, Table 10 shows superiority of Scenario 1. Note that here the cost structure has also changed. In Scenario 1 the blood transfusion depot covers the cost of delivery, whereas in Scenario 2 hospitals cover the cost, such that each hospital covers for its own. Furthermore, Scenario 1 improves distribution of blood products in several ways.(i)Total flow time and transportation costs are minimized, as the blood products distribution problem is modeled as Traveling Repairman Problem where the objective is to minimize the total waiting time of the service to be delivered.(ii)Delays in transportation are minimized. This is attained as special-purpose vehicles are allocated for distributing of blood products, which are considered as emergency vehicles and can overpass traffic limitations.(iii)The highest level of quality of service is maintained compared to Scenario 2. This is because vehicles are highly equipped and staffs are well trained. Obviously, this reduces wrong delivery and improves on time delivery, and most importantly “as ordered” delivery.

According to results (Tables 9 and 10), lower travelling time (equivalently travelling distance) is reached when employing Scenario 1. Thus, it is not very far if we state Scenario 1 outperforms Scenario 2. Finally, considering both technical difficulties, which are much more important than costs (as the problem concerns human lives), and also costs, Scenario 1 is superior to Scenario 2. This is indeed very important, because on a real scale problem, like city of Tehran, improvements are substantial. (Considering both private and public hospitals, currently Tehran holds around 150 hospitals. Also note that we have solved random instances of size 150 hospitals in reasonable amount of time, see Table 7.)

It is worthy to add that currently authors are negotiating with Iran Ministry of Health to design a distribution network of blood products in Tehran. This not only requires software infrastructure, but also hardware ones, most importantly, well-equipped vehicles and well-trained staffs. Here, we believe solving the problem should not be an issue, although the access to the appropriate data is.

6. Conclusion

In this paper, we presented a new application of operations research in healthcare. The problem is to find optimal routes when distributing blood products among hospitals such that the total waiting time of hospitals requiring them are minimized. We studied a real blood distribution problem where 24 hospitals were to be served by a maximum of three vehicles. The proposed solution approaches are a mixed-integer programming formulation and a hybrid heuristic. Computational experiments showed the efficiency of the exact approach in finding optimal solutions for the real case. For the purpose of providing high-quality solutions for much larger instances, several random instances containing up to 150 hospitals and 10 vehicles were generated. The performance of the developed hybrid heuristic on these instances is promising, and near-optimal solutions were reported in a short time (almost 3 minutes on average for these problems). The proposed formulation and solution approach provide a basis for many other home healthcare problems and applications. Currently, the authors are working on these problems and applications.

Appendix

For more details see Table 11.