Safety Times for Multistage Assembly System

1University of Belgrade, Faculty of Mechanical Engineering, Kraljice Marije 16, 11000 Belgrade, Serbia 2Vienna University of Technology, Institute of Computer Languages, Argentinierstr. 8, 1040 Vienna, Austria 3University of Belgrade, Faculty of Forestry, Kneza Viseslava 1, 11000 Belgrade, Serbia 4University of Novi Sad, Faculty of Agriculture, Trg Dositeja Obradovica 8, 21000 Novi Sad, Serbia 5Department of Electrical Engineering, International Islamic University, Islamabad, Pakistan


Introduction
The newsboy problem is well-known.Different variants and extensions of the newsboy problem are used in many applicative scenarios (e.g., [1][2][3][4][5][6][7][8][9]).New requirements, imposed by each specific applicative scenario, introduce more complexity.Therefore, there is a necessity of extended and refined methods in treating of the newsboy problem.Basically, the solution of the newsboy problem is known when it is considered on a simple model.However, it is still an open issue in case of complex models (e.g., if penalties are random values).
In maintenance systems, one of the main problems is uncertainty of events that drastically affect the maintenance efficiency and performance of the system.Typically, these events are random and correlated with demanding time or lead time for some reserve part, or time needed to prepare some reserve part.It is often supposed that a lead time or time needed to prepare some reserve part equals 0 (a reserve part is directly available) or constant, which is a very rough assumption.As many factors influence the above mentioned time, it has a random deviation, which significantly affects the system performance.The classical newsboy problem should be adapted to the given scenario and therefore, the appropriate mathematical models should be developed.
The following assumptions are used: (1) a multistage assembly network is represented by a tree in the graph theory [10], (2) planning of interrelated assembly operations with independent stochastic operation times, and (3) delivery of items in a requested time (taking into account that it cannot be done until all end items are ready).As a consequence, either delay costs appear (if the given delivery requirement cannot be satisfied), or holding costs appear (if end items are ready earlier than the requested delivery time).Therefore, the goal is to find optimal safety times in order to minimize the total expected costs.The considered problem includes several stages, and therefore finding the exact solution is too difficult.In such situations, the usual way of coping with the problem is the application of a certain heuristic that is based on successive applications of the solution of simpler one or two-stage assembly network [11][12][13].A hierarchical structure is approximated by a series of two-level structures.Under certain assumptions, the two-stage problem can be solved as the classical newsboy problem [14][15][16].In [17], the authors consider the problem of contract assembly with combined supply lead time and demand quantity uncertainty and use stochastic operation times as well as stochastic demand.A more general class of stochastic assembly problems is considered in [18].In [12], the authors calculate lower and upper limits for the decision variables and lower and upper bounds for the objective function.In [13], the focus is on dynamic, continuous-time generalization of the singleperiod newsboy problem.They reduce a multilevel system into a series of discrete newsboy problems and give an appropriate algorithm.
Reference [19] solves the newsboy problem by considering multiple discounts and formulates and solves the newsboy problem by considering multiple upgrades as well as a mixed newsboy problem characterized with multiple discounts and upgrades.Reference [20] investigates a multiperiod production planning problem for multiple products, where the uncertain demand and the accumulated effort of market investment on demand are considered.Reference [21] studies a newsvendor model with discrete demand and shows that the optimal ordering decision with discrete demand is very different from that with continuous demand.Reference [22] considers the inverse newsvendor problem as a variant of the traditional newsvendor problem and provides an analysis of the problem under the assumptions of normally and exponentially distributed service times.
The research work presented in [11] inspired us and it serves as our starting point.Namely, on the basis of the mathematical model in [11], we made a generalization.So, [11] is a special case of the presented model (see Section 3, (3)).Our model is built up on the extended assumption that the summed holding costs of individual stage components are greater than or equal to holding costs of the assembled component.The motivation for such a generalization came from the practical use-case (already mentioned maintenance systems).However, there are many similar real-life scenarios (e.g., holding costs of a completely assembled computer that waits to the installation are less than the sum of holding costs of each individual component).Therefore, the introduced relaxation in the starting assumption is important in order to encompass and properly address the abovementioned class of real-life use-cases.
In this paper, we suggest a generalized approach with a random time assumed for the components in leafs (e.g., a delivery from some external source or a preparation for installation).The approach is based on a methodology of decomposition of the general problem into series of simple ones.
The novelty and the main contribution of this paper can be summarized as (i) A new model for determining of critical times in one stage assembly network is given.(ii) A theoretical result given in Theorem 1 (Section 3) is new: the existence and uniqueness of the solution in a given model is proved.
(iii) A heuristic procedure described for the determination of the so-called starting times for all assembly/production steps is new.
(iv) This generalized model can be used in multistage assembly networks.

Problem Statement
In order to better clarify the problem itself, we separately provide the assumptions and notations used in the following text.

Notations.
The following notations are used to establish mathematical models: (ii) In case of several end components, the delay cost is calculated as the maximum delay; (iii) Delay costs associated with other nodes are not taken in consideration, although they may affect the delay at node 0.
Consider node k, and suppose that t k ≥ max (t l +L l , t m +L m , t n +L n ), i.e., component k assembling, cannot start before preceding stage components are not assembled.Holding costs at node k are obtained until the next stage assembly operation starts or the final delivery.The holding costs during the operations are disregarded.
A graphical presentation of the considered assembly network can be seen on Figure 1.It consists of the arcs, which represent the assembling operations, and nodes.Node i is the node where an assembling operation for component i starts, whereas the assembling times L i are independent random variables with continuous cumulative distributions F i and density f i , and e i is positive holding cost per time unit.
The goal is to find optimal safety times, i.e., to plan assembling operations so that the total expected holding and delay costs are minimized [11].

Two-Stage System
As previously mentioned, the approach is based on a methodology of decomposition of the general problem into series of simple ones.An approximate solution is obtained by successively applying the exact solution for a two-stage system.Therefore, we first consider and analyze the system presented in Figure 2.
Let  0 be the optimal starting time for assembly operation at final node 0. Consider stage that precedes the assembly operation at final node 0. The stage delay costs can be expressed as b( −  0 ) + , where b is a positive stage delay cost per time unit,  = max 1≤≤ {  } and T i = t i +L i .Two types of holding costs are incurred.Individual holding costs for component i are incurred until all stage assembly operations are finished and can be expressed as ℎ  =   ( −   ).Stage holding costs are incurred after all stage assembly operations are finished, and can be expressed as ĥ = ê( 0 − ) + , where ê is a positive stage holding cost per time unit.The costs C=C(t 1 ,t 2 ,. ..,t n ) can be expressed as: It is natural (if a delay in the previous stage (with regard to found optimal time t 0 ) exists, from linearity and additivity of the expectation, the expected costs/delay will be transferred to the final stage.) to assume that delay ( −  0 ) + will produce the same delay in the next (final) stage, i.e., corresponding costs should be equal b( −  0 ) + = ( −  0 ) + thus b =  and in the same way ê = .The cost function can be simplified as follows: Notice that the single-stage system (several end items) can be easily obtained as a special case of the two-stage system with final node assembly operation time  ≡ 0 and  = ∑  =1   (see Figure 3).
Especially if n=1,  =  1 , we get the standard Newsboy problem  =  1 (  − 1 ) + +( 1 −  ) + with the optimal solution given by:  1 =   −  −1 1 (/( +  1 )).Consider the general case again.As x + =x+(-x) + , the cost function can be expressed as We assume the following: the summed holding costs per time unit of first stage components are greater than or equal to the holding cost per time unit of second stage component at final node 0 (Figure 2), i.e., ∑  =1   ≥ .For ∑  =1   =  the cost function is which is the model analysed in [11].As   =   +   , i.e., T T = max  {  }, we get If we denote   =  0 −   , the cost function is The expected cost function It can be concluded that the expected cost function Z has a unique minimum in a finite point.
Proof.First, we show that the expected cost function Z given by ( 6) is a convex function with arguments  1 ,  2 , . . .,   .The second part It can be concluded that the exact solution is found in the two-stage system taking into account the introduced assumption.

Optimal Starting Times
thus the expected cost function is: The stationary point of function Z is the solution of the following system of equations: Finally, the stationary point is the solution of the following system of nonlinear equations: The obtained system of equations can be solved by numerical integration.Let us suppose that ( * 1 ,  * 2 , . . .,  *  ) is the solution of the above system.The optimal starting times  *  =  0 −  *  for the two-stage system can be determined accordingly.

An Extension to Multistage System.
Let us consider the two-stage system with root in node i that is a part of the multistage system (Figure 4).The cost function for the stage that precedes the assembly operation at node i is (see ( 2)) where (T i ) j =(t i ) j +(L i ) j is the availability of component (i) j ,   = max  {(  )  } and b ,ê  the corresponding stage delay and stage holding costs.It is necessary to estimate stage delay costs b , i.e., stage holding costs ê .As mentioned in the problem formulation, delays at other nodes may affect the delay at final node 0.
One possible methodology could be to distribute delay costs at final node 0 down the assembly tree.
Assume that the assembly operation at node i starts with delay Δ  .A change of delay costs Δ   at final node 0 is So, the expectation of random variable Δ  /Δ  can be taken as stage delay costs b at node i: where    −  is the distribution function of random variable   −   .On the other side, if stage assembly operations are finished before requested time, that will affect holding costs at final node 0.
Assume that the assembly operation at node i started in some previous moment   −Δ  .A change of holding costs Δ   at final node 0, is The expectation of random variable Δ  /Δ  can be taken as stage holding costs ê at node i: The cost function C i at stage that precedes assembly operation at node i is completely defined: and the expected cost function   = (  ) is as follows: In order to provide the convergence of the above process, it is necessary to show that   is a convex function; i.e., it is enough to prove ∑  (  )  ≥ ê .This condition is valid as it follows directly from the general hypothesis: the summed holding costs per time unit of stage components are greater than or equal to the holding cost per time unit of stage component that follows (Figure 4), i.e., ∑  (  )  ≥   : In this way, the problem solution of safety time estimation can be extended to the multistage assembly network.At the end of this section, we provide a brief sum-up of the proposed procedure.
After finding optimal times in the two-stage tree and a root in final node 0 (12), the summed holding costs and delay costs from final node 0 are being distributed to nodes of the previous stage by using the described procedure (( 15), ( 17)).With each distribution to m nodes of the previous stage (Figure 4), the starting problem is being decomposed into m problems of the lower orders, which are further solved in the same way.

Some Numerical Results and Discussion
We used two sets of sample problems in order to evaluate the suggested approximate procedure.

Problem Sets.
The following problem set is used to evaluate the suggested approximate procedure.The threestage networks are considered in Figures 5 and 6.As it can be seen from the pictures, the problem set 1 is decomposed into two 2-stage problems, whereas the problem set 2 is more complex and decomposed into three 2-stage problems.
We use the following parameters for a numerical example.Holding costs e i at every node are constant and equal to 1, whereas delay costs b at final node 0 belong to {5, 25, 50}.For assembly operations time exponential distribution  : E(1), i.e., normal distribution  : N(0, 1)is used.Further, we assume that distributions of operation times are: operation time equal 1. Dispersion of given distributions is parameterized in the following way.Instead of distribution L, a + (1-a)⋅L, a∈ {0, 0.5, 0.8} is taken, this way, the expected operation time remains 1, but dispersion will be changed.The goal is to determine optimal starting times t 2 , t 3 , t 4 when a certain assembly operations should be activated in order to minimize the total expected costs.Assume that the requested delivery time is t d = 0.

Parameter Sensitivity and Comparison
Analysis.Tables 1  and 2, respectively, show obtained results for 9 test examples.
Table 1 presents the results for the example in Figure 5, whereas Table 2 presents the results for the example in Figure 6.For every test example, there are two subrows.The first subrow contains the values of parameters a and b and the exact solutions for optimal starting times t 2 , t 3 , t 4, i.e. t 3 , t 4 , t 5, and t 6 in which minimum of the total expected costs is reached.The exact solutions are obtained in MatLab.
The second subrow contains solutions obtained by using the proposed approximation and the value of the total expected costs, i.e., deviation of minimal total expected costs expressed in percentage.The column "%" in Tables 1 and 2 presents a relative error between expected costs of the exact solution and expected costs of the solution obtained by using the proposed heuristic.For fixed value of parameter b and with an increased value of parameter a, a relative error decreases.From Table 1, it can be seen that the suggested model treats changes in the assembly system in an adequate way: the smaller dispersion (i.e., the bigger value of parameter a), the less deviation from the exact solution.From the other side, the deviation increases with a disproportion of costs, h and b (h is always equal to 1, b takes values 5, 25, 50).In Table 2, the slightly increased value of error in some cases of the problem set 2 compared to the problem set 1 can be influenced by * ) the increased complexity of the considered problem and * ) the existence of more critical paths L6L2L0, L4L1L0 (as the time distributions L6 and L4 are doubled).
Comparing our results with results of the problem set 2 from Axsater [11] (which is more realistic than the problem set 1), the following could be noticed: (i) regarding deviation of the assembly operation duration both techniques show similar behavior, i.e., a greater deviation of assembly operations durations (smaller value of the parameter a) produce a greater relative error; (ii) in our case, the moments of initiating certain operations on the nodes (i) j preceding the operation corresponding to node i are mainly estimated in same way-either all starting times are overestimated or all of them are underestimated (see Tables 1 and 2), i.e., t 3opt < t 3approx , t 4opt < t 4approx , Table 1.rows (1) , ( 4) , ( 7) , ( 5) , ( 6) , ( 8) , ( 9) Table 2. rows (3) , ( 4) , ( 5) , ( 6) , ( 7) , ( 8) , ( 9) This is expected since our proposed heuristic consists of series of consecutive goals and each goal represents an "overall" optimization of starting times corresponding to nodes () j so that an assembly operation at node  that follows could be started at presumed time   .This is not the case for Axsater [11] procedure, i.e., t s1opt < t s1approx and t s2opt > t s2approx (Table 3).
So, we obtained an improvement of the results presented in [11].The obtained results are promising: (1) the values of deviation are relatively small, and (2) the time improvement becomes significant with the increase of levels; i.e., the greater number of level is introduced in the problem, and the better improvement in results is visible.One of the motivations for the study was also to consider problems with more complex assembly graphs.In a more complex model, it could be noticed a confirmation in the effectiveness of the used heuristic, i.e., a confirmation of the promising results obtained in the simpler model.A comparison of our suggested heuristic with exact solutions on even more complex graphs is not presented as large and complex instances are difficult to be optimally solved by using exact methods [11].
As already mentioned, [11] elaborated that the considered problem with several stages is, in general, too difficult to be solved exactly.So, the important research question that we considered is to construct a heuristic that should be improved by means of two issues: (1) the quality of approximation and (2) a generalization, so that it can "cover" different practical use-cases, which was originally our motivation to start with.In the presented numerical results, we estimated its quality of approximation.The benefit of using the suggested heuristics is the possibility to solve the starting problem for an assembly network with arbitrary many stages.

Statistical Analysis.
In order to draw the conclusion with statistical confidence, i.e., to establish whether the performance on the two selected problem sets differs with statistical significance, the ANOVA test is performed.The comparative analysis showed that, for fixed value of parameter b and with an increased value of parameter a, a relative error decreases.Now, we take another point of view in which we would like to analyze the following: (i) Is the average value of the relative error ("%") significant, for different, variable values of parameters a and b?
(ii) As a consequence of the previous point, what could we expect about the relative error ("%") without regard to specific values of a and b, if problem set changes?
The analysis is set up as follows: (1) Statement of heuristic performance hypothesis: For the two problem sets, the relative error "%" distribution is compared using the null-hypothesis as H 0 : There is no difference in the mean "%" between the two problem sets.H 1 : There is a difference in mean "%" between the two problem sets.
(2) Significance level: The test is performed at a significance level of  = 0.05.
(3) Interpretation of results: If the p-value from the ANOVA-test is less than 0.05, the results are referred to as significantly different and the null-hypothesis, H 0 , is rejected at the defined significance level.If the p-value is greater than 0.05, then the null-hypothesis is not rejected.
According to the obtained data from Tables 1 and 2, the pvalue is 0,3291, which supports the null-hypothesis.

Conclusion
New or changed requirements in maintenance systems increase overall complexity and impose the necessity of finding the solution of appropriate complex models.In such a case, the solution can be found only by using some approximate techniques.This paper considers a multistage assembly system in which the goal is to find optimal starting assembly times taking into account randomness of various events.We proposed one approximation of the multilevel hierarchical structure in the assembly procedure by a series of two-stage structures.A new model for determining optimal starting times in two-stage assembly network is given.The existence and uniqueness of the solution in a given model is proven.Further, this generalized model can be used in arbitrary-level assembly systems.Future work will consider finding the problem solution of safety time estimation in the two-stage assembly network with further relaxation of the introduced assumption (e.g., through an appropriate parameterization, e.g., holding costs and delay costs could be also time functions, random variables).
) is a linear function and thus a convex function, whereas the convexity of the first part and the third part, respectively, (  ∑ =1   − )  (max  {  −   }) ( + )  ((max  {  −   }) and linearity of expectation and convexity of a maximum of convex functions.That means that the expected cost function Z is a convex function as the sum of convex functions.On the other side, as lim ‖‖→+∞  = +∞ (see (1)) from monotonicity of expectation it follows that lim ‖‖→+∞  = lim ‖‖→+∞ (C) = +∞.Properties of the expected cost function Z listed above provide the existence of a unique minimum in finite point ( * 1 ,  * 2 , . . .,  *  ).

Figure 4 :
Figure 4: A two-stage system as a part of the multistage system.

Figure 5 :
Figure 5: An assembly network 1 for the considered problem set.

Figure 6 :
Figure 6: An assembly network 2 for the considered problem set.

Table 1 :
Results for the example in Figure5.

Table 2 :
Results for the example in Figure6.