Performance Analysis of a Kitting Process as a Paired Queue

Nowadays, customers request more variation in a company’s product assortment leading to an increased amount of parts moving around on the shop floor. To cope with this tendency, a kitting process can be implemented. Kitting is the operation of collecting the necessary parts for a given end product in a specific container, called a kit, prior to arriving at an assembly unit. As kitting performance is critical to the overall cost and performance of the manufacturing system, this paper analyses a two-part kitting process as a Markovian model. In particular, kitting is studied as a paired queue, thereby accounting for stochastic part arrivals, and kit assembly times. Using sparse matrix techniques, we assess the impact of kitting interruptions, bursty part arrivals and phase-type distributed kit assembly times on the behaviour of the part buffers. Finally, a cost-profit analysis of kitting processes is conducted and an approximation for a two-part kitting process is established.


Introduction
Nowadays manufacturing systems are often composed of multiple in-house fabrication units [1].The semifinished products stemming from these units are the input materials for other fabrication units or for assembly lines.Hence, efficient transport of materials between the different stages of the production process is key for overall production cost minimisation.Kitting is a particular strategy for supplying materials to an assembly line.Instead of delivering in containers, each holding a single type part and all holding the same number of parts, kitting collects the necessary set of parts for an individual end product in a specific container, referred to as kit, prior to arriving at an assembly unit [1][2][3][4][5][6].
Kitting mitigates storage space requirements at the assembly station since no part inventories need to be kept there.Moreover, parts are placed in proper positions in the container such that assembly time reductions can be realised.Additional benefits include reduced learning time of the workers at the assembly stations and increased quality of the product.Although kitting is a nonvalue adding activity, its application can reduce the overall materials handling time [6].Indeed activities such as selecting and gripping parts are performed more efficiently.Furthermore, the whole operator walking time is drastically reduced or even eliminated since kits, each containing a complete set of components, are brought to the assembly station [7].The advantages mentioned above do not come for free since the kitting operation itself incurs additional costs such as the time and effort for planning the allocation of the parts into kits and the kit preparation itself.Moreover, the introduction of a kitting operation in a production process involves a major investment and the effect on efficiency is uncertain.Therefore, it is important to analyse the performance of kitting in a production environment prior to its actual introduction.This is the subject of the present paper.
In the literature, most authors consider a kitting process as a queueing system with stochastic part arrivals and kit assembly times.Hopp and Simon [8] developed a model for a kitting process with exponentially distributed processing times for kits and Poisson arrivals.They found accurate bounds for the required buffer capacity of kitting processes with two parts.Explicitly accounting for finite buffer capacities, Som et al. [3] further refined the results of Hopp and Simon.
Of course real buffers always have a finite capacity, the capacity being constrained by the storage room.However, if the capacity is large enough, we can have a good approximation of a process with a finite capacity on the basis of a model with unlimited capacity.This means that there is always enough space for upcoming parts, which simplifies the analysis.Unfortunately, the assumption of an infinite buffer is not valid for kitting processes.If the capacity is assumed to be infinite, then the model will degrade to an unstable stochastic system.Harrison [9] showed for a multiple input generalisation of the GI/G/1 queue that it is necessary to impose a restriction on the size of the buffer to ensure stability in the operations of a kitting process.Under this assumption, the probability to have a certain long-term stock position is equal and independent of the current stock position.This was also demonstrated by Latouche [10] who studied waiting lines with paired customers.We can consider this analysis as an abstraction of a kitting process with two types of parts.
In this work, we focus on a kitting process modulated by a Markovian environment.The introduction of this environment allows us to study kitting under more realistic stochastic assumptions: kitting interruptions, bursty part arrivals, phase-type distributed kit assembly times, and so forth.Our paper extends the results on kitting in a Markovian environment [11].
The remainder of this paper is organised as follows.Section 2 describes the kitting process at hand.In Section 3, Chapman-Kolmogorov equations are derived and their numerical solution is discussed.To illustrate our approach, Section 4 considers a number of numerical examples.In particular, we assess the impact of kitting interruptions, bursty part arrivals, and phase-type distributed kit assembly times on the behaviour of the part buffers.Then, a costprofit analysis of kitting processes is conducted and an approximation for a two-part kitting process is established.Finally, conclusions are drawn in Section 5.

Model Description
In this paper, we study a two-queue kitting process, as depicted in Figure 1.Each queue has a finite capacity-let  ℓ denote the capacity of buffer ℓ, ℓ = 1, 2-and models the inventory of parts of a single type.New parts arrive at the buffers and if both buffers are nonempty, a kit is assembled by collecting a part from each buffer.Hence, departures from the buffers are synchronised, and the buffers are paired.Operation of part buffers therefore considerably differs from other queueing systems.
Arrivals at both buffers are modelled by a Markovian arrival process and kit assembly is not instantaneous.For ease of modelling, it is assumed that there is a modulating Markov process, arrival and service rates depending on the state of this process.To be more precise, the kitting process is modelled as a continuous-time Markov process with state space C 1 × C 2 × K, whereby C ℓ = {0, . . .,  ℓ } for ℓ = 1, 2 and with K = {1, 2, . . ., } being the state space of the modulating process.At any time, the state of the kitting process is described by the triplet [, , ], with  and  being the number of parts in the first and second buffer, respectively, and  being the state of the modulating process.We now describe the state transitions.
(i) The state of the modulating process can change when there are neither arrivals nor departures.Let  (iii) Analogously, the state of the modulating process may remain the same or may change when there is a departure (in each buffer).Let   and M denote the corresponding transition rate and generator matrix, respectively.
Summarising, arrivals at and departures from the buffers are described by the generator matrices A, L 1 , L 2 , and M.So far, no diagonal elements of A have been defined.To simplify notation, it will be further assumed that the diagonal elements are chosen such that the row sums of A + L 1 + L 2 + M are zero.The computational method employed here does not require any homogeneity of the generator matrices.When required by the applications at hand, intensities may depend on the buffer content.In this case, we introduce superscripts to make this dependence explicit.For example, M (,)  denotes the generator matrix of state transitions with departure when there are  parts in buffer 1 and  parts in buffer 2. In addition, we use arguments for rates as we already used superscripts and subscripts to distinguish the arrival rates at the different queues.For example,  (ℓ)   (, ) denotes the arrival rate at buffer ℓ = 1, 2 from state  to state  when there are  parts in buffer 1 and  parts in buffer 2.
Example 1.In the most basic setting, parts arrive at the buffers in accordance to an independent Poisson processes with rate  1 and  2 and kit assembly times are exponentially distributed with parameter .In this case, there is no need to have a modulating Markov process, and the state is completely described by the number of parts in each buffer, (, ).We have Example 2. To account for burstiness in the arrival process of the parts at the different buffers, the modulating process allows the mitigation of the Poissonian arrival assumptions.We can replace the Poisson processes by a two-class Markovian arrival processes.Multiclass Markovian arrival processes allow for intricate correlation and can be efficiently characterised from trace data [12,13].As we have two types of arrivals, the Markovian arrival process is described by the generator matrix Λ 1 of transitions with arrivals at buffer 1, the generator matrix Λ 2 with arrivals at buffer 2, and the generator matrix Λ 0 without arrivals.As usual, the diagonal elements of Λ 0 are negative and ensure that the row sums of Λ 0 + Λ 1 + Λ 2 are zero.Retaining exponentially distributed kit assembly times, we have Here I denotes the identity matrix.
Example 3. As for the arrival processes, the model at hand is sufficiently flexible to include phase-type distributed kit assembly times.The phase-type distribution is completely characterised by an initial probability vector  and the matrix T which corresponds to non-absorbing transitions [14].Let t  = −Te  be the column vector with the rates to the absorbing state and let f be a row vector with zero elements except the first one.Assuming Poisson arrivals in both buffers (with rates  1 and  2 , resp.), we get the following matrices: Here, it is assumed that the background state equals 1 if one of the buffers is empty.When service starts again, the background state is chosen in accordance with the probability vector .

Analysis
Having established the modelling assumptions and settled our notation, we now focus on the analysis of the kitting process.

Balance Equations.
We aim to define a set of equations for the steady-state probability vector for the Markov process [ 1 (),  2 (), ()],  ℓ () being the number of parts in buffer ℓ at time  and () being the state of the background process at time .
] be the steady-state probability to be in state [, , ] and let (, ) be the vector with elements   (, ), for  ∈ K. Figure 2 shows a fragment of the transition rate diagram of the kitting model in state [, , ].As mentioned above, two independent input streams arrive at the buffers with intensity  (ℓ)   and are processed into kits with intensity   .Upon completion of a kit, the content of both buffers is decreased by 1.Note that we only show the transitions whereby the modulating Markov process remains in state .Moreover, possible dependence of the transition rates on the buffer sizes is not indicated.
The following performance measures are of interest.
(i) The mean buffer 1 and 2 content is  1 and  2 , respectively, (ii) The variance of the buffer 1 and 2 is Var  1 and Var  2 , respectively, (iii) The effective load of the system  eff is the amount of time that kitting is ongoing.As kitting is only ongoing when none of the buffers is empty, we have (iv) Let the throughput  be defined as the number of kits departing from the system per time unit.Taking into account all possible states from which we can have a departure, we find (v) The shortage probability  is the probability that one of the buffers is empty as  =  (1) (0) +  (2) (0) −  (0, 0) .
(vi) The loss probability   in buffer  is the probability that an arriving part cannot be stored in buffer ,  = 1, 2.
By noting that the accepted arrival load equals the departure rate, we find where   is the arrival load of part  = {1, 2}.If the arrival load in both buffers is the same, then the loss probability is also the same: If the arrival load is not the same, then the excess load in the most loaded buffer is lost as well.

Methodology: Sparse Matrix Techniques.
Queueing models for kitting processes are rather complicated.Indeed, the modelled kitting process has a multidimensional state space.Even for relative moderate buffer capacity, the multidimensionality leads to huge state spaces; this is the socalled state-space explosion problem.For many queueing systems, infinite-buffer assumptions may mitigate this problem.
Given some buffer system with finite capacity, more efficient numerical routines can be constructed for the corresponding queueing system with infinite capacity.Unfortunately, as mentioned above, the infinite-buffer capacity assumption is not applicable for kitting processes and therefore cannot simplify the analysis.Recall that the infinite-capacity model is always unstable.For all input parameters except trivial ones (no arrivals), some or all of the queues grow unbounded with positive probability.Consequently, the multidimensionality of the state space and the inapplicability of the infinite-buffer assumption yield Markov processes with a finite but very large state space.However, the number of possible state transitions from any specific state is limited.This means that most of the entries in the generator matrix are zero; the matrix is sparse.As illustrated by the numerical examples, using sparse matrices and their associated specialised algorithms result in manageable memory consumption and processing times, compared to standard algorithms.
The method used here to solve the sparse matrix equation is the projection method GMRES (generalized minimum residual).Details can be found in Stewart [15] (p.194-197); Philippe et al. [16], Buchholz [17], and Saad and Schultz [18].To solve the linear equation  = , the GMRES algorithm approximates  by the vector   ∈   in a Krylov subspace   = span{, , . . .,  −1 } which minimises the norm of the residual ‖  − ‖.Since the residual norm is minimised at every step of the method, it is clear that it is nonincreasing.However, the work and storage requirement per iteration increases linearly with the iteration count.Hence, the cost of  iterations grows by ( 2 ) which is a major drawback of GMRES.This limitation is usually overcome by restarting the algorithm.After a chosen number of iterations , the accumulated data is cleared and the intermediate results are used as the initial data for the next  iterations [18].
To ensure fast convergence, it is also key to properly choose the initial vector that is passed on to the algorithm.We rely on MATLAB's build in GMRES algorithm and assume a uniform initial vector if no additional information about the solution is known.However, this is often the case.Calculation speed can be further improved as, in practice, performance measures are not calculated for an isolated set of parameters.For example, when a plot is created, a parameter is varied over a range of values.In this case, a previously calculated steadystate vector for some set of parameters can be used as a first estimate of the steady-state vector for a new "perturbed" set of parameters.Using previously calculated steady-state vectors is trivial if the state spaces corresponding to the parameter sets are equal.In this case, the previously calculated steadystate vector can be passed on unmodified.If the state space changes, the steady-state vector must be rescaled to the new state space.In general, adding zero-probability states if the state space increases or removing states if the state space decreases turns out to be ineffective.This is easily explained by a simple example.Assume that we increase the queue capacity of one of the part buffers.Typically, even for moderate load, a considerable amount of probability mass can be found for queue size equal to capacity.Increasing the queue size and assigning zero probability to the new states is not a good estimate for the new steady-state vector.Also for the system with higher capacity, a considerable amount of probability mass can be found when the queue size equals the capacity (while zero probabilities were assigned).

Numerical Results
With the balance equations at hand, we now illustrate our numerical approach by means of some examples.

Bursty Part Arrivals.
As a first example, we quantify the impact of production inefficiency on the performance of a kitting process.To this end, we compare part buffers with Poisson arrivals to corresponding kitting systems with interrupted Poisson arrivals.The arrival interruptions account for inefficiency in the production process.Kit assembly times are assumed to be exponentially distributed with service rate equal to one, this value being independent of the number of parts in the different buffers.This is a kitting process with Markovian arrivals as described in Example 2 in Section 2.
The interrupted Poisson process considered here is a two-state Markovian process.In the active state, new parts arrive in accordance with a Poisson process with rate  whereas no new parts arrive in the inactive state.Let  and  denote the rate from the active to the inactive state and vice versa, respectively.We then use the following parameters to characterise the interrupted Poisson process: Note that  is the fraction of time that the interrupted Poisson process is active, the absolute time parameter  is the average duration of an active and an inactive period, and  is the arrival load of the parts.Figure 3 shows the mean number of parts in buffer 1 versus the arrival load, for various values of the buffer capacities  1 and  2 for Poisson arrivals (for both buffers) as well as for interrupted Poisson arrivals (again for both buffers).The arrival load is set to  = 0.8 for all curves.In addition, we set  = 0.4 and  = 10 for the interrupted Poisson processes.Clearly, the mean buffer content increases as the arrival load increases as expected.Moreover, if more buffer capacity is available, it will also be used: the mean buffer content increases for increasing values of  1 =  2 .Comparing interrupted Poisson and Poisson processes, burstiness in the production process has a negative impact on performance-more buffering is required-if the queues are not fully loaded ( < 1).As for ordinary queues, the opposite can be observed for overloaded buffers.In this case, the effect of a large burst is mainly reflected in loss and not in additional queue content.Burstiness also yields larger periods without arrivals during which the buffer size decreases.
By numerical examples, we can quantify expected buffer behaviour-for example, more production yields higher buffer content, higher buffer capacity mitigates the loss probability, and so forth.However, less trivial behaviour can be observed as well.Figure 4 depicts the probability that Mean number of parts in buffer 1 the buffer is full versus the buffer capacity  1 =  2 .We compare performance of kitting with Poisson arrivals to kitting with interrupted Poisson arrivals at one buffer and at both buffers.As in the preceding figure, the interrupted Poisson processes are characterised by  = 0.4 and  = 10.
As expected, the probability that the buffer is full decreases for increasing values of the buffer capacities.Moreover, to reduce this probability, more buffer capacity is required for the case of two interrupted Poisson processes than for the case of two Poisson processes.For the kitting process with one Poisson and one interrupted Poisson process, non-trivial performance results can be observed.Namely, interruptions in the production of a part more negatively affect buffer performance of the other part.Indeed, buffer 1 is full with higher probability if the arrivals at buffer 2 are interrupted than if the arrivals at buffer 1 are interrupted.First note that the loss probabilities in both buffers are the same.For the Poisson buffer, parts are rejected at the arrival rate  1 whenever the buffer is full.For the IPP buffer, parts are rejected when the buffer is full and the arrival process is in the on-state at rate  2 /.As the loss probability in both buffers is the same, we have Pr As the second queue is more likely to be full when there are arrivals we have Pr

Phase-Type
Distributed Kit Assembly Times.The second numerical example quantifies the impact of the distribution of the kit assembly times on kitting performance.In particular, we here study Erlang-distributed kit assembly times.
Limiting ourselves to Poisson arrivals to both buffers, this numerical example fits Example 3 of Section 2.   Figures 5 and 6 depict the mean number of parts in buffer 1 and the loss probability in buffer 1 for the kitting process and, as a reference point, for the M/E/1/N queue as well.In both figures, the arrival load is varied and different values of the variance of the kitting time distribution are assumed as indicated.The mean kitting time is equal to 1 for all curves and the capacity of both buffers is equal to 20.In underload ( < 1), kitting performs worse than the M/E/1/N queue: the mean buffer content and the loss probability have a higher value.This follows from the fact that kitting stops when one of the buffers is empty.By increasing the load, it is obvious  that the buffer content converges to the capacity and the loss probability to one.It is most interesting to observe that the shape of the service time distribution only has a small effect on these performance measures.Indeed, there is no significant performance difference when  2 equals 1/4 and when it equals 1/8.

Cost and Profit Analysis.
We now add a cost structure to the kitting process under study.In particular, cost and profit for the kitting systems of Sections 4.1 and 4.2 are analysed.
The proposed cost function is where  1 is the holding cost of a part in the buffer,  2 is the shortage cost in one or in both of the buffers, and  3 is the loss cost.Note that for all figures, the input parameters are symmetric for both parts such that  1 =  2 and  1 =  2 .

Poisson Arrivals.
In Figure 7, the total cost of applying kitting systems versus the buffer capacity (varying from 1 to 30) is depicted.Limiting ourselves to Poisson arrivals for both buffers,   =   = 0.8 for part  = {1, 2}.We compare kitting performance for different costs as depicted.As expected, higher cost values lead to higher total cost.The cost structure also affects the optimal buffer capacity.When looking at the different costs separately, we observe that a higher holding cost decreases the optimal buffer capacity (from 15 to 12) and a higher loss cost increases the optimal capacity (from 12 to 14), whereas the optimal buffer capacity remains the same (12) for a higher shortage cost.Obviously, buffering is interesting when the storage cost is low and when the cost of rejecting parts (because of a full buffer) is high.It is most interesting to observe that, as the capacity increases, the cost models converge when the sum of the holding and the shortage cost ( 1 +  2 ) is equal.Indeed, as the loss tends to zero, the loss cost tends to zero as well.Furthermore, when the capacity equals one, the state space of the kitting model has size 4 and hence the cost is easily calculated explicitly as Figure 9: Shape of the service time distribution does not impact significantly the value of the optimal buffer capacity.
by a sale unit equal to 100.Assuming a maximum storage room of 30 parts, we observe that the optimal buffer capacity is 5, 11, and 12 for the depicted arrival processes, respectively.Consequently, kitting systems under production inefficiency require much higher storage space, especially when profit is applied as the parameter to determine the optimal buffer capacity.Moreover, the optimal capacity is very sensitive to the burstiness parameters  and .In the example at hand, the plots suggest that the cost function is a convex function of the buffer capacity.However, one can also easily choose the cost parameters to obtain nonconvex cost functions.shows that GMRES performs considerably is better than LU decomposition for capacity.

Conclusion
In this paper, we evaluate the performance of two-part kitting processes in a Markovian setting.Furthermore, a costprofit analysis is conducted and an approximation for our model is given.Note that the particularity of the studied kitting systems is that the part buffers are paired.This means that each demand requires both parts such that a kit can only be assembled if both inventories are nonempty.Methodologically, we have applied sparse matrix techniques (e.g., GMRES) as most of the entries in the generator matrix have a value equal to zero.The solution is not exact but performs well in terms of solution speed and accuracy.
As our numerical results show, the interplay between the different queues leads to complex performance behaviour.For example, interruptions in the production of a part more negatively affect buffer performance of the other part.Overall, we observe extreme sensitivity of kitting performance on arrival process parameters while performance is reasonable insensitive to variation of the kitting time distribution.Furthermore, we determine the optimal buffer capacity based on a cost-profit analysis and establish a quite good approximation for our two-part kitting process.

2 MFigure 1 :
Figure 1: Kitting process: the buffers are on the left, and the triangularly shaped kitting process is on the right.

Figure 3 :
Figure3: When the queues are not fully loaded, burstiness in the production process has a negative impact on the performance.

Figure 4 :
Figure4: Irregularity in the production of a part leads to a higher probability to have a full buffer for the other part.

Figure 5 :
Figure 5: Variation of the kitting time distribution does not have a significant impact on the mean number of parts.

Figure 6 :
Figure 6: Variation of the kitting time distribution does not have a significant impact on the loss probability.

Figure 7 :Figure 8 :
Figure 7: Each type of cost has a different impact on the value of the optimal buffer capacity.
equal to 2, a shortage cost  2 equal to 55, and a loss cost  3 equal to 40.On the left figure, the optimal buffer capacity for Poisson arrivals is 12, for an interrupted Poisson arrival at one buffer it equals 22, and for an interrupted Poisson process at both buffers the optimal capacity is 28.As expected, the optimal buffer capacity increases as the burstiness in the production process increases.Concerning the profit analysis, we assume that the yield equals the throughput multiplied ) 4.3.2.Bursty Part Arrivals.Next, we analyse the cost and the profit of kitting systems with different Markovian arrivals.We