Determining Optimal Replacement Policy with an

Full text

Translate

Turn on search term navigation

[ProQuest: [...] denotes non US-ASCII text; see PDF]

Shengliang Zong 1 and Guorong Chai 1 and Yana Su 2

Academic Editor:Marco Mussetta

1, School of Management, Lanzhou University, Lanzhou, Gansu 730000, China

2, College of Economics and Management, Lanzhou Institute of Technology, Lanzhou 730050, China

Received 22 April 2016; Revised 23 November 2016; Accepted 15 December 2016; 18 January 2017

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Power equipment such as a wind power generation system is operating in a dynamic and uncertain environment. Such equipment may fail due to random events such as sand/rain storms, ice rain, heavy snow, and extreme temperatures. The service interruptions caused by the equipment failure will affect the company's performance in terms of revenue and customer satisfaction. For example, during January 10 to February 5, 2008, an unexpected strong East Asian winter monsoon occurred with extremely low temperatures, blizzard conditions, and freezing rain in southern China. Due to its unusual persistence and intensity, this extreme event caused extensive damage and general disruption over southern China. Most of the infrastructure damages were attributed to freezing rain, which led to a huge amount of broken power equipment, electric power transmission lines, and chaotic traffic conditions for southern China. The disastrous heavy snow caused direct economic losses of 53.8 billion RMB [1, 2].

In practice, the power equipment may be interrupted by random shocks, which represent severe weather conditions such as thunder storm and freezing rain or operational contingencies such as short circuit. When a shock arrives, it will cause a random amount of damage to the system. When the accumulated amount of shock damage exceeds a specified threshold, the system will fail. Therefore, it is necessary to consider the maintenance problem under random shocks. Chen et al. [3] investigated the impact of the voltage imbalance on the reliability of power distribution systems. Kiureghian [4] considered the problem of distribution of electric substation equipment suffering from earthquake.

Recently, shock models were utilized to model the operating time. Conceptually, a system fails due to the shock effect on the system. While most of existing shock models were based on the accumulated or extreme damage causing a system failure. Li [5] first introduced the δ -shock model to avoid measuring amount of damage which may not be easy in many situations. Although the shock's damage amount is difficult to measure, its arrival time can be accurately recorded. Therefore, the δ -shock model focuses mainly on the frequency of shocks rather than their accumulated damage amount of shocks. In a δ -shock model, a shock causes a system failure which also called a deadly shock if the time interval between two successive shocks is smaller than a critical value δ . The threshold δ is usually a constant. In this paper, we adopt a general δ -shock model by letting δ be an exponentially distributed random variable with parameter varying with number of repairs [6]. The power equipment after repair can be more fragile and more prone to fail again. As a result, the threshold of a deadly shock is increasing with n , which is the number of repairs.

There are extensive studies on the maintenance problems for the system with operating and repair stages. This is mainly because some classical assumptions are not realistic in modeling the real systems. These assumptions include repaired system becoming "as good as new" and a dialed system being replaced by a new one immediately. Barlow and Proschan [7] introduced an imperfect repair model, where the repair is prefect with probability p and minimal with probability 1-p. We also consider the nonzero repair and replacement times in contrast to a classical assumption that repairs are instantaneous. In most practical situations, to reflect the aging process of the system, the consecutive repair times are assumed to become longer and longer till the system is replaced with a new one according to some replacement rule. The geometric processes (GP) method, with its practical advantages in characterizing the whole operation process of the components in a system [8], gains worldwide scholar's attention in recent years, and many researches have been carried out based on GP in different application domains. Lam [9] first introduced GP to study the maintenance problem for the deteriorating system. Finkelstein [10] generalized Lam's work by scale transformation. Most of the studies based on GP denote a good applicability in the related fields, but their assumption that the property of the component after premaintenance is as good as new is not fully compatible with the electrical power equipment power [11].

Ensuring system function and life are two very important maintenance goals [12]. Li et al. [13] designed an age-related dynamical model to create a timely maintenance schedule that extends the life of the equipment and improves its efficiency. Preventive maintenance also needs to consider how to reduce costs [14]. Castanier et al. [15] proposed a parametric maintenance decision framework to coordinate inspection and replacement of the two components and to minimize the long-run maintenance cost of the system. Lin and Huang [16] used a nonhomogeneous Poisson process with a power law failure intensity to describe the deterioration of a repairable system and to minimize the expected total cost per unit time to determine the optimal nonperiodic maintenance schedule. As we all know, the goal of the maintenance decision is not only to minimize costs. We adopt the replacement policy by satisfying certain requirements of equipment availability. Farag and Al-Baiyat [17] consider two major problems of power system design optimization. One is to maximize the reliability under the cost constraint and the other is to minimize the cost under the reliability constraint. Levitin and Lisnianski [18] studied the optimization problem of the power system containing a bridge substructure in which redundant elements are included in order to provide a desired level of reliability. There are more research works on power system structure optimization under reliability constraints, for example, Lisnianski et al. [19] and Bevrani et al. [20].

Existing research goal is to minimize the total cost or the long-run maintenance cost of the system, and the objective function mainly is linear. In this paper, the objective function is an integer nonlinear programming, and constraints are an integer nonlinear inequality. In order to try to solve this kind of problem in an appropriate way, a method based on genetic algorithms (GA) is proposed to explore the optimality of this problem. GA is heuristic search techniques, which are full of promise to deal with complex optimization problems [21]. GA can also solve optimization problems possessing arbitrary degrees of nonlinearity, discontinuity, and stochasticity. Levitin and Lisnianski [22] used GA to minimize the sum of costs of system modernization actions over the study period while satisfying reliability constraints at each stage. Abdulwhab et al. [23] used GA optimization technique to maximize the overall system reliability for a specified future time period in which a number of generating units are to be removed from service for preventive maintenance.

In this paper, we propose a δ -shock model, which will be helpful for practitioners to find the appropriate machine maintenance and replacement policies, for electrical power equipment. We assume that the shocks will arrive according to a renewal process with an exponential distributed interval time of shocks and that the threshold of a deadly shock forms a geometrical increase. Furthermore, the distribution function of the threshold follows exponential distribution. The assumptions and analysis of the model are introduced in Section 2. Section 3 presents an optimization model which is to minimize the long-run average cost per unit time subject to some availability constraints. The GA is introduced in Section 4. Finally, Section 5 provides a case study using our approach.

2. Model and Analysis

2.1. Model Assumptions

To construct the δ -shock model, we first introduce the definition of GP as follows.

Definition 1.

A stochastic process {_ξn , n=1,2,...} is called a geometric increasing (or decreasing) process, if there is a real number 0<a<=1 (or a≥1), such that {^an-1_ξn , n=1,2,...} forms a renewal process. a is called the ratio of the GP.

Moreover, let E(_ξ1 )=τ and Var(_ξ1 )=^σ2 . Then, we have [figure omitted; refer to PDF]

Therefore, a, τ, and ^σ2 are three parameters for a GP.

We then propose a repair model for electrical power equipment under availability constraints based on the following assumptions.

(I) At the beginning, a new system is installed. Whenever the system fails, it will be immediately either repaired or replaced by an identical system.

(II) The system in operation is subject to environment shocks which arrive according to a Poisson process with rate _λ1 or E_Xn =1/_λ1 . Let _Xn (n=1,2,...) be the time interval between the (n-1)th and the nth shocks. Then _Xn is assumed to be exponentially distributed with the distribution function denoted by F(x). Let _δn be another exponentially distributed random variable associated with _Xn . We assume that the sequence {_δn , n=1,2,...} forms an increasing geometric process with 0<a<=1. Then, _δn has cumulative distribution function Q(^an-1 x), where Q(x) is the cumulative distribution function of _δ1 which is a brand new system's threshold of the shock model. Assume that E_δ1 =1/_λ2 . {_Xn ,_δn } follows a δ-shock model if the system fails at nth shock which satisfies _Xn <=_δn , and then the life time or equivalently the operating time is the sum of all _Xn until the one satisfying the above condition. Further, we assume that _Xn is independent of _δn .

(III) Let _Tn be the operating time between the (n-1)th and nth repairs. {_Tn , n=1,2,...} is a sequence of stochastically decreasing random variables.

(IV) Let _Yn denote the repair time after the nth failure. Then the distribution function of _Yn is G(^bn-1 y), where 0<b<=1 is a constant. Let E(_Yn ) denote the expected repair time after the nth failure. Assume that E(Y)=μ>0.

(V) _Tn and _Yn , n=1,2,..., are independent.

(VI) The reward and cost structure is assumed as follows: the reward rate when the system is operating is r. The repair cost rate of the system is c. The replacement cost consists of a fixed replacement cost R and variable cost v=_rp Z, where _rp is the variable cost rate and Z is the replacement time with a mean of E(Z)=t.

(VII) The N-replacement policy is adopted. Under such a policy, the system will be replaced with a new and identical one at the Nth failure.

(VIII) The equipment must meet the minimum availability requirements. Let p represent the availability of the equipment. In other words, the proportion of the nonproductive time, for example, repair or replacement, should not be too high.

2.2. Maintenance Cost

Since _Yn , n=1,2,..., forms an increasing GP with ratio 0<b<=1, we have [figure omitted; refer to PDF]

Thus, the maintenance cost function is shown as [figure omitted; refer to PDF]

2.3. Operating Revenue

Let E(_Tn ) be the expected operating time of the system between the (n-1)th and nth repairs, and let _lni be the interval time between the (i-1)th and ith shock following the (n-1)th repair. Then the number of shocks experienced until the next failure can be written as [figure omitted; refer to PDF]

We have [figure omitted; refer to PDF]

Here _Mn follows a geometric distribution. That is [figure omitted; refer to PDF] where _pn =P (a shock following the (n-1)th repair is a deadly shock) and _qn =1-_pn .

Therefore, we have E_Mn =1/_pn . As _Mn is a stopping time with respect to the random sequence {_lni , i=1,2,...}, which are independent identically distributed random variables. According to Wald's equation [24], we have [figure omitted; refer to PDF]

Based on Assumption (II), and the fact that F(x) and Q(x) are exponential distributions, we have [figure omitted; refer to PDF]

Furthermore, as _lni and _δn are independent and have the marginal exponential distributions with means of 1/_λ1 and 1/^an-1_λ2 , therefore, we obtain [figure omitted; refer to PDF]

Then [figure omitted; refer to PDF]

Consequently, [figure omitted; refer to PDF]

Then, the operating revenue of the system, denoted by OR, is given by [figure omitted; refer to PDF]

2.4. Replacement Cost

Based on Assumption (VI), the replacement cost function, denoted by RC, can be written as [figure omitted; refer to PDF]

2.5. Long-Run Average Cost

Let C(x) denote the long-run average cost per unit time of the system. According to renewal reward theorem [24] [figure omitted; refer to PDF]

Let W be the expected length of a renewal cycle of the system under replacement policy N. Thus, we have [figure omitted; refer to PDF]

Let TC be the total expected cost of the renewal cycle of the system under replacement policy N. Then, [figure omitted; refer to PDF]

Consequently, [figure omitted; refer to PDF]

3. Optimization

3.1. System Availability

Using (13) and (17), the availability function of the system can be written as [figure omitted; refer to PDF]

3.2. Minimizing the Average Cost

The optimal policy can be found by minimizing the cost of each unit time subject to a minimum level of system availability. The optimization problem is [figure omitted; refer to PDF]

Thus, we can see that the optimization problem is an integer nonlinear programming problem. The analytical solution is not effective to solve this optimization problem, because the objective function C(x) is a nonconvex and extremely complex function. Several classical computational techniques, for example, branch and bound technique, cutting planes technique, implicit enumeration, and out approximation, which are reasonably efficient, have been proposed in literature for solving integer nonlinear programming problems [25-27]. These techniques are applicable to a particular class of problem. In the case of nonconvex problems these techniques may cut-off the global optimal.

4. Genetic Algorithms

In the last two decades, a number of stochastic algorithms have been proposed for solving the integer nonlinear programming problems. Genetic algorithms have been proved to be a useful tool in solving global optimization problems [28]. GA is general purpose population based stochastic search techniques which mimic the principles of natural selection and genetics. It can be implemented easily with the computer programming. In particular, GA is more effective to solve the large scale, real life, discrete, and continuous optimization problems without making unrealistic assumptions and approximations. Recently, GA has been proposed to solve integer programming problems by Soltani and Shafiei [29] and Wong et al. [30]. Applying GA to solve nonlinear mixed integer optimization problems was reported by Pal et al. [31] and Deep et al. [32].

Our objective function is an integer nonlinear programming and constraints are an integer nonlinear inequality. For solving such a problem, we have developed a real coded GA for integer variables.

It is generally accepted that any GA consists of the following basic components.

Step 1 (set the parameters).

The parameters of the GA in our problem include population size, probability of crossover (pc), probability of mutation (pm), and stopping conditions. The population size should be neither too large nor too small. We set the parameters as follows:

Population size = 10, pc = 0.90, pm = 0.10, stall generations as the stopping conditions.

Step 2 (generate initial populations).

An initial population of size n is randomly generated. Any general procedure to get them can be applied. Nevertheless, in the implementation, we prefer to build the initial population with feasible chromosomes. This device guarantees that, upon termination, the method will provide a bilevel feasible solution. A positive integer is randomly generated according to the nature of the problem. It follows the next law l=ceil (min x (1) + (max x (1) - min x (1)) [low *] rand (Size, 1)).

Step 3 (calculate the fitness values).

Each chromosome in the population is evaluated by creating the estimated value for each coefficient in (19). The fitness value of each chromosome is then obtained from (20). Let, f(i) be the fitness value for each chromosome and the total fitness value of the populations is T, where T = sum[...] f(i).

Step 4 (perform selection).

The roulette wheel selection mechanism is the selection strategy chosen in our paper. Each time, a single chromosome is selected from the current generation to create a new generation. The selection process is as follows. Each time, a random number h is generated from the range [0,1]. Let p_r(i) be the weight of the ith chromosome in the whole populations, where p_r(i) =f(i)/sum f(i). If h<p_r(i) , then this chromosome is selected. Otherwise, the jth chromosome is selected if p_r(j-1) <h<p_r(j) (i<j). This selection process is continued until the new population has been created. The chromosome which has the higher fitness value tends to be picked out to the next generation.

Step 5 (perform crossover).

Crossover is the key to the power of genetic algorithms. The purpose of the crossover is to generate offspring which endorses them with the common character of their parents. Calvete et al. [33] proposed two types of crossover methods. One is called "variable-to-variable" crossover and the other is called "basis-to-basis" crossover. Lin [34] used the one-point crossover method, which randomly selects one cut-point and exchanges the right parts of two randomly selected parents in the population to generate offspring. The crossover method of this paper uses also the one-point method.

Step 6 (perform mutation).

Mutation is an operator that acts on various chromosomes by randomly exchanging one of its indices. Mutation resets a selected position in a chromosome to a randomly generated real number in [0,1]. The mutation rate is usually a very small value; on average, 0.3% of the total populations undergo mutation [34]. The mutation method used here is the midvalue method, where the mutation of chromosome took place in the middle of the population. The new chromosome contains the indices of the new basis.

Step 7 (return to Step 3 until the stopping condition is met).

There are many stopping conditions for GA, including generations, time limit, fitness limit, stall generations, and stall time. The programs will stop until one of the stopping conditions is met. The stopping condition used here is stall generations.

5. Empirical Results

In this section we analyze an example using the proposed GA. All experiments were performed in MATLAB on a PC with an Intel Core i7 2.80 GHz and 4 G RAM. The data is simulation data which are used to further evaluate the proposed approaches. The simulation data are generated from the analytical curves. Table 1 shows the value of the parameters.

Table 1: Value of the parameters.

Parameters	a	b	c	µ	R	_{r p}	t	_{λ 1}	_{λ 2}	p
Value	0.95	0.94	2	4	8000	3	20	0.005	0.02	0.98

Figure 1 shows the relations between C(x) and x. Note that the black dot meets constraints and the red dots do not. It is obvious that the best results can be obtained at x=32.

Figure 1: The plot of C(x) against x.

[figure omitted; refer to PDF]

In order to simplify the calculations, we set the lower bound to be 0 and the upper bound to be 1000, which is a positive infinity in our problem. The appropriate parameters of GA here are set as follows.

: Population size = 10

: Probability of crossover = 0.90

: Probability of mutation = 0.10

: Stall generations = 200

The stopping condition for the GA is to reach the maximum stall generations. The simulation result is shown in Figure 2. It is interesting to observe that, in initial generations, the GA finds the configuration to minimize the cost. The fitness is suddenly reduced from -0.2 to -2.2 on the 10th generations. After 10 generations, the fitness has decreased. The global fitness is obtained at about the 45th generation. Best fitness is -2.4950 and the corresponding variable is 32.

Figure 2: Best fitness terms evolution during the GA generations (t=1.6719 s).

[figure omitted; refer to PDF]

In addition, the GA is repeatedly executed 10 times to ensure that the global fitness can be achieved so that the local optimization can be avoided. As is shown in Table 2, it is observed that we can find the running times of different executions are different, and the 5th stop and the 7th stop do not reach the global optimal results.

Table 2: Result of 10 times (p=0.98).

	Value	Individual	Time (seconds)	Y / N
1st stop	-2.4950	32	1.6250	Y
2nd stop	-2.4950	32	1.6406	Y
3rd stop	-2.4950	32	1.6406	Y
4th stop	-2.4950	32	1.6406	Y
5th stop	-2.4872	30	1.6563	N
6th stop	-2.4950	32	1.6719	Y
7th stop	-2.4872	30	1.6875	N
8th stop	-2.4950	32	1.6719	Y
9th stop	-2.4950	32	1.6875	Y
10th stop	-2.4950	32	1.7031	Y

Furthermore, we simulate the fitness with different constraints. Figure 3 shows that the three lines have the same trend of converging to the optimal fitness in the 10th generation. Table 3 shows that the number of repairs is decreasing and the long-run average cost per unit time is increasing in p. As is shown in Table 3, when the availability of the equipment is 0.97, the best fitness is -2.5026 and the corresponding variable is 37. Meanwhile, if the availability of the equipment is increasing to 0.99, the best fitness will become -2.3261 and the corresponding variable will become 17. In practice, the availability of the equipment is decreasing with the repair times. Thus, you have to maintain more times until the equipment reached its availability, and you will obtain more benefits due to the equipment operating longer times.

Table 3: Comparison of solutions with different p values.

	p = 0.97	p = 0.98	p = 0.99
Individual	Fitness	Individual	Fitness	Individual	Fitness
Value	37	-2.5026	32	-2.4950	17	-2.3261
5	-0.3%	0	0%	-15	6.8%

Figure 3: Best fitness terms evolution with differences p during the GA generations.

[figure omitted; refer to PDF]

6. Conclusions

In practice, the power equipment may be interrupted by random shocks, which represent severe weather conditions such as thunder storm and freezing rain or operational contingencies such as short circuit. When the accumulated amount of shock damage exceeds a specified threshold, the system will fail. Therefore, it is necessary to consider the maintenance problem under random shocks. In this paper, we present a nonlinear integer optimization model for power equipment to minimize the expected long-run average cost rate under availability constraints. Such a model is useful for system managers to find the optimal maintenance policies. We proposed to use the GA to solve the optimization problem. GA appears to be a strong method for solving complex global optimization problems in which gradient related methods usually give local optima. GA has global scope over the search space so they get global optima. Such an approach can be very effective in solving this type of optimization problems with respect to cost and availability. As we all know, GA is an approximate algorithm. The optimal solution may be different due to it is approximate optimal solution each time. And it may take more times until we obtain the optimal replacement policy when the scale is larger. Meanwhile, the results show that the number of repairs is decreasing and the long-run average cost per unit time is increasing with availability. In other words, you have to maintain more times before reaching the availability of the equipment, and you will obtain more benefits due to the equipment operating longer times. There is not only the availability but also the maintenance cost and effects on the optimal replacement policy. In a word, our model provides a useful quantitative tool for managers to evaluate the system performance and design an optimal maintenance policy.

In this paper, the shock model process is assumed to follow a Poisson process. Extending such a shock process into a more general nonhomogeneous Poisson process or renewal process can be a good topic for future research. Meanwhile, we only considered the average cost minimization problem. A possible future work also can be to maximize the availability under the cost constraints. In addition, our approach can be extended to study multiequipment of the parallel or unparallel systems.

Notations