1 INTRODUCTION
Emergency service systems (ESSs) still see workforce scheduling as a challenge and have even been studied within operations research/management science for more than fifty years (^{Green & Kolesar, 2004}; ^{Simpson & Hancock, 2009}). The main objective of workforce scheduling is to minimize operational costs while maintaining adequate service level. However, the daily operations of ESSs are surrounded by randomness and require studies to consider variations in the morning, afternoon, and night operations (^{Schmid, 2012}; ^{Souza et al., 2015}).
Due to the nature of the events they answer, ESSs usually must have a low workload through the day. This is a characteristic of systems with a low event frequency (^{Chiyoshi et al., 2011}). The literature is thin for such systems. The majority of the literature studies high event frequency systems such as call centers.
The majority of the queueing models used on employee scheduling are stationary or piecewise constant (^{Green et al., 2001}). This approximation creates distortions related to the user arrival process and, more rarely, to the service process (^{Kim & Whitt, 2014}). There are also distortions to model the process of servers leaving the system at the end of their shifts. (^{Ingolfsson et al., 2007}).
The use of timedependent queueing models makes employee scheduling a larger but more realistic challenge by revealing its nonlinear and dynamic behavior. These characteristics are associated with problems involving difficulttofind optimal solutions, so the use of heuristics becomes necessary, as raising server requirements is often unrealistic (^{Ingolfsson et al., 2010}). We can enumerate some efforts on how the genetic algorithm is used (^{Ingolfsson et al., 2002}), as an integer programming heuristic (^{Ingolfsson et al., 2010}), and as the heuristic used in ^{Chung and Min (2014}).
Shiftscheduling problems still face limitations. As mentioned before, the use of stationary approximations is extremely common (^{Defraeye and Van Nieuwenhuyse, 2016}). Also, ^{Ingolfsson et al. (2002}, ^{2010}) pointed out that improved heuristics with lower computational efforts are desired, as the use of genetic algorithm is computationally expensive. Also, ^{Feldman et al. (2008}) and ^{Defraeye and Van Nieuwenhuyse (2016}) stated that systems with low arrival rates are a challenge for many optimization methods, turning smallscale systems into avenues for further research. Finally, ^{Defraeye and Van Nieuwenhuyse (2016}) indicated that the development of easytocompute problems could benefit the study of more complex performance metrics and their relations in queueing systems.
Hence, this paper aims to improve the shiftscheduling process revisiting the police precinct study from ^{Kolesar et al. (1975}) and ^{Ingolfsson et al. (2002}). Firstly, the genetic algorithm used in ^{Ingolfsson et al. (2002}) is analyzed regarding its convergence, fitness function, number of iterations, and computational times. The purpose of the genetic algorithm here is to determine how many servers will work on each possible shift. This answer must satisfy all constraints and have adequate cost and delay probability. This application of the genetic algorithm for shiftscheduling problems is based on the application found in ^{Ingolfsson et al. (2002}). We also propose a new fitness function and initial population generation method based on previous results, to make the genetic algorithm more reliable. A local search heuristic is also proposed for the same purpose. Such a simple heuristic is passed through an initial test, and its results are compared to that of the genetic algorithm through a sensitivity analysis decreasing and increasing arrival rates in a smallscale system. Therefore, this paper contributes to the literature by improving the genetic algorithm for shift scheduling problems, by presenting a simple and fast local search heuristic, and by evaluating these methods in a benchmark case study.
The remaining paper is structured as follows. Section 2 presents the related literature to the topic. Section 3 presents the methodology of our study, including: data used for the New York police precinct; the queueing model and its performance measures used to evaluate schedules; the base genetic algorithm and the proposed improvements; and the proposed local search heuristic. Section 4 shows the results of the proposed methods in the benchmark case and a brief analysis about both the genetic algorithm and proposed local search heuristic under demand changes. Finally, Section 5 shows some final considerations about this study and perspectives for further research.
2 RELATED LITERATURE
^{Dantzig (1954}) and ^{Edie (1954}) mark the advent of the personal scheduling problems back in the 1950s. Firstly, ^{Edie (1954}) collected data on tolls to obtain hourly personnel requirements and suggested that a linear programming problem could make use of such requirements to calculate an optimal shiftscheduling. The calculations for the requirements used stationary queueing models as the M/M/s (Erlang C) for each day hour. ^{Green et al. (2001}) called this approach as Stationary Independent PeriodbyPeriod (SIPP). ^{Dantzig (1954}) presented an integer programming solution for ^{Edie’s (1954}) suggestion. The combined formulation of ^{Dantzig (1954}) and ^{Edie (1954}) are known as SIPPIP (SIPP solved with Integer Programing) method (^{Ingolfsson et al., 2002}).
^{Baker (1976}) classified personal scheduling in three categories: shiftscheduling, daysoff scheduling, and tourscheduling (combination of the first two categories). Our focus is on shiftscheduling problem, where planning considers a daily planning horizon for the scheduling problem (^{Van den Bergh et al., 2013}).
The literature on shiftscheduling evolved a lot since its advent in the 1950s. A quick search in a Scopus database revealed 121 published papers about shiftscheduling in the last 20 years. Half of them were published after 2012. Preferred journals include (but are not limited to) the European Journal of Operational Research, Journal of Scheduling, Computers and Operations Research, and Annals of Operations Research.
Building a shift schedule has some steps (^{Ingolfsson et al. 2002}). One can find papers on the literature with different orders and method to get the schedules (see ^{Thompson (1997}) and the references therein). Typically, it is necessary to (^{Buffa et al., 1976}):
Forecast periodbyperiod demand,
Convert demand forecasts into staffing requirements,
Determine set of possible shifts,
Find a shift schedule that obey all requirements raised and minimize a goal as number of employees or labor costs.
Demand forecasts are made using historical data (^{Gillard and Knight 2014}; ^{Angelo et al., 2017}); converting calls (arrivals) into arrival rates obeying a nonhomogeneous Poisson process (^{Ingolfsson et al. 2002}). ^{Kim and Whitt (2014}) shows statistical tests to evaluate if an arrival rate is wellmodeled using nonhomogeneous Poisson processes.
The requirements often respect a service level measure, which we consider as the probability that there is going to be at least one idle server when a new user arrives in the system. ^{Gillard and Knight (2014}) uses Singular Spectrum Analysis (a technique to decompose a time series into oscillatory components and noise) to generate demand forecasts and requirements because they consider that other means do not capture seasonal variations properly. ^{Defraeye and Van Nieuwenhuyse (2016}) and ^{Schwarz et al. (2016}) presented comprehensive literature reviews on this topic.
Possible shifts depend on the problem modeled. For instance, typically, industrial companies have nonoverlapping shifts (^{Van den Bergh et al., 2013}). Most of the problems are related to systems with a great number of events, such as call centers (^{Brown et al., 2005}; ^{Dietz, 2011}; ^{Kim & Ha, 2012}). Some other papers are related to hospital emergency units (^{Creemers et al., 2012}; ^{Green & Kolesar, 1991}; ^{Green, et al., 2001}; ^{Gillard & Knight, 2014}; ^{Sungur et al., 2017}). One can find other services in ^{Ko and Gautam (2010}), ^{Chung and Min (2014}), Thompson (1993, ^{1997}), and ^{Jennings et al. (1996}). Among papers about ESS, one can cite ^{Kolesar et al. (1975}), using data from a police precinct.
As the problems themselves, finding a solution (schedule) that obeys all requirements and optimizes the goal is a great challenge and many methods arise in the literature. ^{Alfares (2004}) categorized solutions of scheduling problems in 10 groups: manual solutions, integer programing, implicit modeling, decomposition, goal programming, working set generation, linearprograming based solution, construction/improvement, metaheuristics, and others.
^{Ingolfsson et al. (2002}) developed a method for workforce scheduling. The approach modeled the structure of possible shifts and the dynamic nature of demand. A genetic algorithm solved a benchmark case (shown in Section 3.1) and a numerical solution of nonstationary queueing system evaluated service levels. The authors found that the SIPPIP method can significantly overestimate the service level that results from a schedule.
The nonstationary queueing model used in ^{Ingolfsson et al. (2002}) is a set of the ChapmanKolmogorov equations (^{Kleinrock, 1975}). However, ^{Ingolfsson (2005}) and ^{Ingolfsson et al. (2007}) pointed flaws on the direct use of ChapmanKolmogorov equations without considering what happens the moment a busy server is scheduled to leave. The literature identified two possibilities for this moment, which are called endofshift disciplines. Firstly, the service is preempted and the user rejoins the head of the queue (preemptive discipline). Secondly, the server completes the job before leaving (exhaustive discipline). ^{Ingolfsson (2005}) points that the exhaustive discipline is often more realistic when users are humans. However, the exhaustive discipline still requires more advances to represent mealtime breaks (^{Beojone, 2017}). Those interested on nonstationary queueing models may find detailed explanation for them and their numerical solutions in ^{Schwarz et al. (2016}), ^{Defraeye and Van Nieuwenhuyse (2016}) and the references therein.
3 METHODOLOGY
This section presents the studied case, a benchmark problem from the literature. It also presents the schedule evaluator (nonstationary queueing model). Further in this section we present the genetic algorithm and the proposed improvements for a shiftscheduling problem. Finally, we present the proposed local search heuristic.
A portion of our study is inspired on ^{Ingolfsson et al. (2002}). Therefore, it is imperative to clarify the points of convergence between this paper and ^{Ingolfsson et al. (2002}). Firstly, we used the same benchmark problem from New York’s police. Secondly, we used the same queueing model. Finally, we use a similar construction for the genetic algorithm, which we improve in the next subsections. Besides the improvements on the genetic algorithm, this paper also contributes with a new local search heuristic, compared to ^{Ingolfsson et al. (2002}).
3.1 Data and studied case
^{Kolesar et al. (1975}) studied a police precinct in New York to develop a shiftscheduling problem based on stationary queueing models. ^{Ingolfsson et al. (2002}) revisited this case to implement a shiftschedule problem using a genetic algorithm, which is explained in more detail ahead. Table 1 summarizes the main information about the New York police case seen in ^{Kolesar et al. (1975}) that will be used in our study, such as the arrival rates and the working hours for each possible tour. Note that there are three major shifts, with four mealtime breaks between 2 and 5 hours after the shift’s beginning, for a total of 12 possible shifts, as used in ^{Ingolfsson et al. (2002}). Service time is considered constant throughout the day, with a constant rate µ equal to 2 calls/hour. This information is the same as used in both ^{Kolesar et al. (1975}) and ^{Ingolfsson et al. (2002}).
t  0  1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23 

λ (t)  9.8  9.6  8.7  7.6  6.7  5.3  4.1  3.2  2.5  2.5  2.9  3.8  4.3  5.0  5.9  6.6  7.8  8.6  9.4  9.8  10.2  10.4  10.2  10.0 
Tour presence  1  1  0  1  1  1  1  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 
1  1  1  0  1  1  1  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  
1  1  1  1  0  1  1  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  
1  1  1  1  1  0  1  1  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  
0  0  0  0  0  0  0  0  1  1  0  1  1  1  1  1  0  0  0  0  0  0  0  0  
0  0  0  0  0  0  0  0  1  1  1  0  1  1  1  1  0  0  0  0  0  0  0  0  
0  0  0  0  0  0  0  0  1  1  1  1  0  1  1  1  0  0  0  0  0  0  0  0  
0  0  0  0  0  0  0  0  1  1  1  1  1  0  1  1  0  0  0  0  0  0  0  0  
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  1  0  1  1  1  1  1  
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  1  1  0  1  1  1  1  
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  1  1  1  0  1  1  1  
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  1  1  1  1  0  1  1 
3.2 Schedule evaluator (nonstationary queueing model)
According to ^{Gillard and Knight (2014}), it is usual to ESSs to operate under stochastic conditions, so it is absolutely important to be able to plan for such situations. Also, these stochastic conditions are timevarying, and it is reasonable to approximate an ESS with a M(t)/M/s(t) queueing model. For those unfamiliar with Kendall’s notation, M(t)/M/s(t) represents a queueing system with a nonhomogeneous Poisson arrival process M(t), a negative exponential service time M, and a timevarying number of working servers s(t). Arrival and service processes are represented by the λ(t) and µ rates, respectively, in calls/hour.
As servers have defined shifts, it is necessary to define what happens in case a server is providing service the moment he is scheduled to leave. The behavior assumed by servers in this circunstance is called endofshift discipline. As ^{Ingolfsson et al. (2010}) did, we adopted preemptive endofshift discipline to facilitate the comparison of our computational results. Hence the nonstationary queueing model used in this paper is the same as the one seen in ^{Ingolfsson et al. (2002}), Green and Soares (2007), and ^{Ingolfsson et al. (2010}). The remaining of this subsection presents a brief explanation on the queueing model.
Let π _{ n } (t) denote the transient state probability that n customers are in the system at time t. We can write the ChapmanKolmogorov forward differential equations for a timedependent multiserver queue as seen in Equation (1) to describe the exact behavior of such a system. The ChapmanKolmogorov forward differential equations were solved, numerically, using the RungeKutta method.
In this paper, nonstationary queueing models are firstly aimed at calculating instantaneous delay probabilities and hourly mean delay probabilities. To do so, we define delay probability π _{ D } (t) as the probability that an arriving new user must wait before being served. Our objective is to guarantee a lower π _{ D } (t) than a given probability P for all values of time t. Equation (2) shows the delay probability as calculated in nonstationary queueing models.
We also need to calculate the hourly mean service level, defined as follows in Equation (3):
As λ(t) is constant on the interval
We use the schedule evaluator (nonstationary queueing model) to support the decision process of both heuristics from subsections 3.3 (genetic algorithm) and 3.4 (local search).
3.3 Genetic Algorithm
Many characteristics of shiftscheduling problems lead us to the use of heuristics to solve them, such as the size and complexity of the search space and the nonlinear behavior of delay probabilities, along with its discontinuities. We started this application with a standard genetic algorithm.
^{Goldberg (1989}) presented the genetic algorithm. Its main idea is based on Darwin’s evolutionary theory, so it needs a population of individuals that will compete to select the most adapted individual. Each individual might be selected to breed with another one. The most adapted individuals will have higher chances to breed and will thus pass its characteristics to the next generations. The least adapted individuals will have smaller chances to breed and will thus tend to disappear. Individuals are also susceptible to random mutations, altering their chances to survive and breed (for more information, see ^{Goldberg, 1989}).
The algorithm uses chromosomes to represent individuals, and each individual is a potential solution for the problem. The most adapted individual to have passed through the population at any time will be chosen as the final solution. These solutions are encoded in the form of binary vectors, as shown in Figure 1. With a fixed number of individuals (called the population), a new generation is created as an offspring of the previous population (no members of the previous generation are included in the next generation). To select individuals (chromosomes  Figure 1) for breeding, we measure their adaptation by calculating a fitness function. Then, two distinct individuals are selected to breed with replacement. Higher fitness leads to higher chances of being selected. Regarding the breeding process, selected individuals have their chromosomes pass through a process called crossover, in which their chromosomes are combined to create two new chromosomes (individuals). After crossover, these new chromosomes may suffer “mutation” at any point, where a bit may change from 0 to 1 and viceversa. When a new generation emerges, we have an iteration, and the selection process begins once more until a number of iterations is reached or another stop criterion becomes true.
The seed for the initial population is the solution of the scheduling problem using SIPP method and integer programing (^{Dantzig, 1954}). This solution method is referred in the remaining of the paper as “SIPPIP method” and the solution itself (shift schedule for servers) from this method is referred as “SIPPIP solution”.
The process for generating the initial population used here was proposed by ^{Ingolfsson et al. (2002}). The initial population size is 60 individuals. The initial population was generated starting from the SIPPIP solution. The SIPPIP was one of the members of the initial population, and the remaining 59 individuals were generated randomly as follows: (1) we calculated the mean and variance of the SIPPIP solution (number of servers in each shift is the piece of data for this calculation); (2) we sampled from a normal distribution with the calculated mean and variance for each shift (defining new shift schedule for all shifts); and (3) sampled numbers were rounded to the nearest integer between 0 and 7 (limits of the adopted chromosome binary encoding; see Figure 1). The mutation probability for each bit was set to 1%, and the maximum number of iterations was set to 100.
Minimizing labor costs here means fewer vehicles, and maintaining an adequate service level means to keep the delay probability lower than 10%. Thus, our fitness functions were built using the same three parameters as in ^{Ingolfsson et al. (2002}): (1) CARS, or the sum of all vehicles from each tour; (2) MAXPROB, or the maximum delay probability along the time horizon; and (3) NUMHOURS, or the number of hours with a higher delay probability than the desired 10%.
The first fitness function used here is based on the third fitness function presented in ^{Ingolfsson et al. (2002}), and it is presented as Equation (5). M indicates a number large enough that the fitness function will not become negative. This fitness function is referred as FITNESS0 in the remainder of the paper. The second stop criterion (besides the number of iterations) for FITNESS0 was a difference between the fitness of the current population’s best individual and average fitness of the population lower than 0.5.
The fitness values of individuals were linearly scaled to ensure that the average fitness value would not change (^{Ingolfsson et al., 2002}). Hence, the selection of individuals for breeding used the scaled value. The purpose of scaling fitness values was to avoid premature convergence during the first iterations and to avoid a random walk among the mediocre individuals in late iterations (^{Goldberg, 1989}). Scaling process is also useful to avoid some bias caused by using a value of M in the fitness function.
3.3.1 FITNESS function choice
However, FITNESS0 puts too much attention to MAXPROB, which cannot indicate whether the service level constraint is respected or not. In the other hand, NUMHOURS is an interesting option, because it only assumes values different from 0 if the service level is higher than the constraint limit. Therefore, as an option for FITNESS0, we worked on a fitness function that tries to obey the servicelevel constraint at all costs. Equation (6) shows FITNESS1, built using the same parameters as before, except for MAXPROB.
3.3.2 Initial population generation
When evaluating the population along the iterations, the breeding process comes across several individuals that are certainly not optimal. However, due to scaling and the nature of genetic algorithms, the chromosomes of these individuals might pass through generations. In a Markovian process, as the queueing models we use, it is very unlikely that a situation in the beginning of the first shift will interfere in the third shift, and even in the second shift. For this reason, we propose another method to generate the initial population, focused in the shifts that do not obey the service level constraint. Before generating any initial population members, we ran the nonstationary queueing model to find the tours in which the delay probability constraint was not obeyed. For example, if the constraint was not obeyed somewhere between 2 and 3 a.m., then all tours beginning at midnight and ending at 8 a.m. were chosen because police officers are in the middle of their shift (mealtime breaks were considered as part of their shifts). This strategy was used to modify only tours that did not obey the service level constraint. Equation (7) shows the proposed method, where
3.4 Local Search Heuristic
If we consider that the SIPPIP solution is a reasonably good initial guess to the problem, it is straightforward to think about using a local search heuristic to find a solution. The module uses the SIPPIP solution and tries to satisfy the delay probability constraints through minor changes to the current solution (for an introduction to local search algorithms, see ^{Luke, 2013}). The main idea is to try to change servers’ tours to satisfy the delay probability constraint. In case it is not possible to fulfill the constraints only changing servers’ tours, we add one new server to a chosen tour and try to manage it until the delay probability constraint is obeyed. In the other hand, if the SIPPIP solution fulfills the delay probability constraint, it will be replaced by a poorer SIPPIP solution, which is calculated after relaxing the delay probability constraint. Table 2 shows the notation used for the local search heuristic.
Notation  Description 

W  Set of all tours. 
w, y  Tours of work (one of the rows from Table 1). 
k  Current iteration. 

Number of servers scheduled to tour w at iteration k. 

Starting time of the mealtime break of tour w. 
WS _{ w }  Set of tours with the same starting time as tour w. 
rC _{ D }  Relaxed delay probability constraint. 
C _{ D }  Delay probability constraint. 
maxt  Observable time horizon. 

Delay probability on iteration k at epoch t. 

Mean delay probability on iteration k at period

Table 3 illustrates the proposed local search heuristic through a pseudocode. It starts on iteration k ← 0 and calculates a solution from the SIPPIP. At this point, the SIPPIP considers a rC
_{
D
} to ensure that C
_{
D
} is still active on the nonstationary queueing model (measured on the schedule evaluator). Then, the code set maxt as 24 hours. Note that the primary stop condition is to make C
_{
D
} inactive. The code verifies whether maxt is acceptable (capable of identifying a
4 COMPUTATIONAL RESULTS
This section presents the results computed for the aforementioned studied case with the genetic algorithm and the proposed improvements. Ahead in this section we present the results for the local search and some comparisons between both solution methods in a sensitivity analysis. All computations ran in a Windows 10 Home, with a processor Intel Core® i5 7300HQ 2.50GHz and 8GB of RAM memory.
4.1 Genetic algorithm
Figure 2 shows the first results for configuration explained here. To obtain more consistent data, we ran the genetic algorithm and calculated an average value for each measure. Hence, Figure 2 does not show actual values for any run of the algorithm but an average value for all of them. Max., min., and avg. are the highest, lowest, and average FITNESS0 values for the population in each iteration, respectively. In observing the reference values of Figure 2, one can see a slight improvement to the average and minimum values for FITNESS0. However, the highest FITNESS0 value is in the initial population (iteration 0) and then falls and increases again, but it never reaches the initial solution value for FITNESS0. On 25 out of 30 runs, the best solution found was the SIPPIP solution. In four runs, the genetic algorithm found a solution with 28 CARS, and in one run, it found a solution with 27 CARS. None of the 30 solutions obeyed the delay probability constraint. It is important to note that the SIPPIP solution found (number of servers for each shift) was different from that shown in ^{Ingolfsson et al. (2002}) but with the same CARS value. This SIPPIP solution had a MAXPROB of 0.12, against 0.14 for the SIPPIP solution from ^{Ingolfsson et al. (2002}).
Figure 3 shows the average results for the genetic algorithm after using FITNESS1 for 30 runs. Differently from FITNESS0, one can see a greater gap between the reference values in the first iterations and greater improvements to the average and minimum values. However, the upper reference values show a slight increase on the first iterations and an even slighter decrease until the final iterations. Only four solutions were found after the 20 first iterations, and 19 were found until iteration 10.
Figure 4 shows the results for both fitness functions. The circle sizes represent the number of times that a solution was the solution found. Firstly, FITNESS0, as mentioned before, found only three solutions, none of which obeyed the delay probability constraint. Secondly, FITNESS1 found 18 different solutions, but it did not obey the delay probability constraint just twice. The most common solution had 31 CARS and a MAXPROB of 0.09; this solution was found eight times.
It is also important to note that changing the fitness function had an effect on runtimes. While FITNESS0 had an average runtime close to 400 seconds, FITNESS1 had 600second runtimes. There was also a difference in variation, with FITNESS0 having greater variance and FITNESS1 having a smaller one. We stored previous solutions in a database, which was erased at each new run. FITNESS0 found its solutions after looking at 2,500 solutions, while FITNESS1 looked at 3,200 solutions, on average.
Figure 6 shows the average reference results for FITNESS1 after using the alternative method for generating the initial population. Note that the initial population had good reference values, even for the poorest solution. In 20 out of 30 runs, the best solution was found among the initial population and the first iteration. This solution uses 30 CARS and has 0.09 MAXPROB. For comparison, using the previous initial population generation method, this solution was found only eight times, versus 26 times using the new method. Also, no solutions were found after the eighth iteration. The upper solutions decreased their FITNESS1 along the iterations, even though the average solutions had only slight improvements.
4.2 Local search heuristic
Figure 7 presents the MAXPROB and CARS results for each iteration of the proposed local search heuristic under different delay probabilities. Each marker represents the results of a single iteration. Larger markers indicate the final solution, and the numbers linked to them denote at which iteration the final solution was found. Hence, no more than 6 iterations separate the SIPPIP solution (iteration 1) from the final solutions. This solution is final, and the local search will not find any other solution for this problem because there were no random choices or steps. Note that, although the algorithm tries to change servers from tours in the best way (taking a server with a mealtime break during a busy period and put it into a tour with a calm mealtime break), it always found worse solutions than the SIPPIP solution. In the other hand, after adding new servers, rearranging them among tours decreased MAXPROB. Moreover, it is interesting to note that the simple idea of adding servers based on their mealtime breaks is a good strategy to reduce delays.
4.3 Comparisons
This final set of tests was aimed at comparing the results from the genetic algorithm and the proposed local search algorithm. The comparisons focused on increases and decreases in demand. The service times remained the same.
Figure 8 presents the number of vehicles for both the genetic algorithm and the proposed local search heuristic as well as the solution from the SIPPIP, for illustrational purposes. The genetic algorithm is represented by circles. The larger the circles are, the more times the genetic algorithm found that number of vehicles in the final solution. This does not mean that these solutions are the same; rather, it only means that they have the same number of vehicles. One can observe that the local search heuristic usually found the same solution as the genetic algorithm while the SIPPIP overestimated the performance, allocating fewer vehicles.
Figure 9 explores the delay probability calculated for final solutions from the genetic algorithm, the proposed local search heuristic, and the SIPPIP. In none of the cases did the SIPPIP solution have acceptable performance. Also, as the demand increased, the genetic algorithm found more nonacceptable solutions (delay probability higher than 0.10). The genetic algorithm only bested the proposed local search heuristic, for a 60% demand increase, in which the proposed local search needed 43 vehicles to guarantee an adequate delay probability, the genetic algorithm found a solution with only 42 vehicles.
Another important point regards runtimes. Figure 10 shows the correlation between the average runtimes and the number of iterations for the proposed local search heuristic. The proposed local search heuristic had almost constant runtimes, much smaller than those seen for the genetic algorithm (about 600 seconds, as seen in Figure 5).
5 FINAL CONSIDERATIONS
We revisited the police precinct study from ^{Kolesar et al. (1975}) and ^{Ingolfsson et al. (2002}) to improve the shiftscheduling process using nonstationary queueing models to evaluate schedules. The genetic algorithm used in ^{Ingolfsson et al. (2002}) received improvements over the fitness function and initial population generation method. We also proposed a local search heuristic, in the form of a hillclimb heuristic. This simple heuristic passed through an initial test, and its results were compared to those of the genetic algorithm through a sensitivity analysis by decreasing and increasing arrival rates on a smallscale system.
Our work has important implications for managers to improve the service quality and efficiency of an ESS. Firstly, the improved genetic algorithm found solutions that obey the delay probability constraint more often and have close to the optimum number of vehicles. Two important features were responsible for these improvements: the new fitness function, which forced the choice of schedules that obey the delay probability constraint, and the new initial population generation method. We illustrated through a series of tests how these new features affect convergence and the genetic algorithm’s final solutions. In addition, results from the improved genetic algorithm led us to the proposed local search heuristic. This heuristic also finds close to optimum schedules with much lower computational expense, especially under low arrival rates. For comparison, the improved genetic algorithm found solutions after 600 seconds, on average, while the proposed local search heuristic found them in less than 5 seconds (even less than 1 second, in some tests). Both the genetic algorithm and local search heuristic determined the appropriate number of vehicles allocated in each possible tour, while considering 1hourlong mealtime breaks.
Differently from a genetic algorithm, the developed heuristic will always come to the same results, as it does not rely on random choices. Furthermore, it finds one final solution from one initial solution (the SIPPIP solution), rather than from a population of solutions. As the proposed local search heuristic does not use chromosomes, its solutions are not limited to the limits that the chromosome encoding might impose on the problem.
The developed heuristic works with only one welldefined goal, making it simple and straightforward to implement. It also relies on the intuitive idea of managing servers’ tours based on their mealtime breaks. Nevertheless, the code for the heuristic is simple enough to accept changes and cope with multiple objectives.
There is a wide field for further research. Moreover, using ChapmanKolmogorov equations is computationally expensive, so we encourage the use of approximate methods like randomization/uniformization (^{Grassmann, 1977}) and the stationary backlogcarryover approach (^{Stolletz, 2008}). The search for better and more efficient heuristics and algorithms is encouraged. Also, there is still room for improvement to the genetic algorithm, such as by using approximations instead of exact models or better fitness functions using aggregate measures. We suggest the use of the proposed local search heuristic for smallscale systems, when considering the use of complex performance measures, to take advantage of the simplicity of our method.