A Condition-Based Maintenance Policy for a Continuously Deteriorating Multi-Unit System with Aperiodic Inspections

(1)

A Condition-Based Maintenance Policy for a Continuously

Deteriorating Multi-Unit System with Aperiodic Inspections

(2)

Master’s Thesis Econometrics, Operations Research and Actuarial Studies Specialization: Operations Research

University of Groningen Minou C.A. Olde Keizer Groningen, November 29, 2012 Supervisors:

Prof. Dr. R.H. Teunter (RUG) G. de Boer (ORTEC)

Co-assessor:

(3)

A Condition-Based Maintenance Policy for a Continuously

Deteriorating Multi-Unit System with Aperiodic Inspections

Minou C.A. Olde Keizer (1774905) November 29, 2012

Abstract

In condition-based maintenance (CBM) (i.e., a predictive maintenance strategy), it is intended to perform maintenance right before a failure occurs by estimating the pending moment of failure based on monitoring a certain condition, such as vibration or temperature. In this Master’s Thesis we consider the CBM optimization approach of Castanier, Grall and B´erenguer (2005)1. Their approach is advanced, compared to those proposed by others, in that it optimizes the inspection moments as well as the condition thresholds (on which the planning of maintenance actions is based) simultaneously. It considers a discrete time system with two machines that operate in series, where the objective is to minimize the long run average maintenance cost per time unit. Our main contributions are as follows. First, we modify the analysis of Castanier et al. (2005) such that all possibilities of their model are included. Second, we analyze an adapted version of their system where two machines operate in parallel rather than in series, and provide new insights on CBM for systems with redundancy.

1

(4)

(5)

Nomenclature

ζi Restore threshold for component i

ζ₁i Opportunistic repair threshold for component i ζi

2 Opportunistic replacement threshold for component i

ξi_l Inspection threshold for component i, l ∈ {0, . . . , Ni− 1}

ξi_N

i Preventive repair threshold for component i

ξi_N

i+1 Preventive replacement threshold for component i

π(x1, x2, . . . , xn) Stationary law of the maintained system state at the start of an

inspec-tion Ai

∞ Long run system availability for component i

C∞ Long run average operating costs

C(t) Cumulative operating costs up to time t Ci

C(t) Corrective replacement costs for component i up to time t

CF(t) Unavailability cost due to failure up to time t

C_Fi (t) Unavailability cost of component i due to failure up to time t CI(t) Inspection costs up to time t

C_Pi (t) Preventive replacement costs for component i up to time t C_Ri (t) Repair costs for component i up to time t

cc Cost of a corrective replacement (per item)

cf Unavailability cost rate (per item and per unit of time)

ci Cost of an inspection (all items)

cp Cost of a preventive replacement (per item)

cr(x0) Cost of a repair with improvement x0 (per item)

cs1 Set-up cost for a repair

cs2 Set-up cost for a replacement (preventive or corrective)

DF(t) Time spent in failed state during period t

Di_F(t) Time spent in failed state during period t for component i ¯

DF(t) Upper bound of time spend in failed state during t

¯

Di_F(t) Upper bound of time spend in failed state during t for component i fi Probability density function of the deterioration increments of

compo-nent i

f∆y Probability density function of the repair in a component’s deterioration

level

Hi(k|y; l) Probability that the unavailability time of component i is k, given a

previous deterioration level of y, with ξ_{N −l}i ≤ y < ξi N −l+1

Li (Fixed) failure level of component i

M Number of subintervals in extended midpoint rule

Ni Number of inspection thresholds for component i

n Number of components in the system

Ppq(r) Probability that both ¯D1F(S) and ¯D2F(S) equal r given conditions p and

q

P_pqi (r) Probability that ¯D_Fi (S) equals r given conditions p and q S Length of a semi-regeneration cycle in steady state

SO1(t) Savings in set-up costs from opportunistic repairs up to time t

(8)

2 NOMENCLATURE

X_ti Condition of component i at time t

(9)

3

1 Introduction

Technical systems often deal with increasing wear and tear caused by usage, age or random shocks. If ignored, this deterioration may eventually cause a system breakdown, which can lead to high costs, system unavailability and safety hazards. Performing preventive mainte-nance can help to prevent failures and its corresponding detriments by repairing or replacing a component before a system breakdown occurs (Gertsbakh (1977)). An advantage of preven-tive maintenance is the freedom to choose the exact point in time of performing a maintenance action on forehand rather then waiting for a failure to occur. In this way it is possible to use moments of low utilization of the system for performing maintenance activities. Common maintenance actions include

Repair

Improving the current state of the system to a better one, not necessarily as good as new. The extent to which the system is improved may be stochastic.

Replacement

Changing the state of the system to the ‘as good as new’ state.

Note that a repair or a replacement can be performed either preventively or correctively. In case of a preventive maintenance action the exact timing depends on the maintenance strategy chosen, while a corrective maintenance operation is performed after a failure has occurred. The different maintenance actions can be combined to form a maintenance strategy. Over the past decades a lot of research has been performed in the field of maintenance strategies, on which a number of surveys has been written, e.g. McCall (1965) and Wang (2002). In Gits (1992) and Barlow and Hunter (1960) it is explained that several maintenance strategies are available, which are schematically represented in Figure 1.

(10)

4 1 INTRODUCTION

A more detailed explanation of the different maintenance strategies is given below. Failure-based maintenance

Replacing the system in case of a failure without performing any preventive maintenance actions.

Block replacement

Performing maintenance after a prespecified operating time independent of the num-ber of intervening failures. If a failure occurs the system is replaced, but preventive maintenance is not rescheduled.

Age-based maintenance

Performing maintenance after a prespecified period of use. When a failure has occurred the system is replaced and preventive maintenance is rescheduled.

Condition-based maintenance

Initiating a maintenance action upon reaching a prespecified condition threshold, and replacing the system in case of a failure. A measurable condition should be available which gives an good indication of the moment of failure (Geurts (1983)).

(11)

5

model by including aperiodic inspection. The system can be inspected at equidistant times, but does not need to be inspected at each possible moment. An advantage of using aperiodic inspection times is that it takes into account that when a system is relatively new inspections are required less often than when the system is approaching the end of its life. This can reduce the costs substantially. In van Noortwijk (2009) it is stated that also Castanier, B´erenguer, and Grall (2003) consider a model which is similar to that of Abdel-Hameed (1987), in which aperiodic inspections are included. The authors restrict the model to a fixed failure level, and include non-negligible durations of maintenance activities. Also the possibilities of repairs are included, and it is assumed that the improvement of the condition due to repair is stochastic. In order to assess the performances of the model the authors consider two performance crite-ria for the proposed maintenance policy, namely the long run system availability and the long run expected maintenance costs. Castanier, Grall, and B´erenguer (2005) adapt the model of Abdel-Hameed (1987) such that it considers two components. Also for this two-unit series system aperiodic inspections are included, and the model is restricted to a fixed failure level. Due to the increasing complexity of the model that arises from considering two components the authors did not include durations of maintenance activities nor repairs. As a performance criterion only the long run average maintenance cost is considered.

Most existing literature on condition-based maintenance considers either a model that op-timizes the critical level at which a preventive replacement should be initiated given the (periodic) inspection intervals or optimizes the inspection intervals given the critical level. The combination of the two, so optimizing the maintenance policy with respect to both the inspection intervals and the critical preventive maintenance level, has only been studied by Grall et al. (2002) and Castanier et al. (2003, 2005) and by Wang (2000). In addition Grall et al. (2002) and Castanier et al. (2003, 2005) consider aperiodic inspection times, which has also been ill-researched. The authors explain that the performances of their models are promising, even better than the classical preventive maintenance methods, although this has not been tested with real data.

In this Master’s Thesis we focus on the discrete time condition-based maintenance models described by Grall et al. (2002) and Castanier et al. (2003, 2005). In Castanier et al. (2005) a two-unit series system is considered, and one of our goals is to apply this model to a system with two identical units that are functioning parallel instead of in series. In this way the system will not completely fail as soon as at least one of the components fails as is the case in the article. We will assess the performances of our model, and compare these with those of classical preventive maintenance strategies. The current model for a two-unit series system is very basic, it does not allow for repairs of the components. We would hence also like to extend the model to this more realistic case, by making use of Castanier et al. (2003). Furthermore, the options of extending the model such that it includes durations of maintenance activities and such that it includes three components will be considered. While reading and extending on the work published by Castanier et al. (2005) we observed that the authors did not always take into account the possibilities of performing opportunistic replacements. Therefore, we present a modified version of this model as well.

(12)

6 1 INTRODUCTION

(13)

7

2 System Description

Below we first describe our deterioration model, after which we explain the maintenance strategies investigated. Here we give an overview of the different maintenance actions that can be performed and we explain the multi-threshold policy that we apply. Furthermore the different models that we consider and compare are defined in detail.

2.1 Deterioration Model

In this section the stochastic deterioration model that we use is described, which is similar to the one defined in Grall et al. (2002) and Castanier et al. (2003, 2005). The discrete time system under consideration consists of n components (n ≥ 1), which independently suffer from increasing wear. Suppose that the condition of component i at time k can be described by a random variable X_ki, with i ∈ {1, . . . , n} and k ∈ N. A failure of component i will occur as soon as its deterioration level exceeds a preset failure level Li, but can only be

noticed upon inspection. This means that the component remains in the failed state until the next planned inspection, an assumption which applies for example in a system where failure does not imply a system stop, but where the quality of produced items is reduced. Also in case of a standby system this assumption is useful. It is assumed that inspections and maintenance actions can only be performed on the system at the discrete equidistant times kδ. The length of δ is usually determined by taking into account available working hours, fluctuating demand for usage of the component or environmental regulations. We set the time unit to δ, i.e., δ = 1, from which it follows that an inspection can be performed at time k (k ∈ N). Although the moment of failure can be estimated more accurately with smaller time units, this assumption is made to reduce the complexity of the model. Note that the degradation of the global system is given by (Xk)k∈N = (Xk1, . . . , Xkn)k∈N. Since

we assume that the degradation of each component will increase over time, we require the random deterioration increments ∆_(k,k+1)Xi _{in each time interval to be nonnegative. This}

(14)

8 2 SYSTEM DESCRIPTION

Figure 2: Example evolution of the system deterioration level of component i.

Let fi denote the probability density function of the deterioration increments ∆(k,k+1)Xi of

component i, for i = 1, . . . , n and for all k ∈ N. From the assumption that the deteriora-tion increments are stadeteriora-tionary and exchangeable it follows that the distribudeteriora-tion funcdeteriora-tions fi

should be infinitely divisible (Lawless and Crowder (2004)). This is for example the case for all gamma distributions. The property of infinitely divisible increments makes a gamma pro-cess suitable for describing deterioration caused by continuous use (Singpurwalla and Wilson (1998)). To be more precise, it is very suitable for describing the steady evolution of wear between the start-up period and the wear accumulation at the end of the system’s life. The choice of a gamma process is confirmed by van Noortwijk (2009), in which it is stated that it is the most appropriate choice for our deterioration model.

2.2 Maintenance Strategy

In this Master’s Thesis we consider the condition-based maintenance models as described in Grall et al. (2002) and Castanier et al. (2003, 2005) and a number of extensions of these models. A summary of the properties of each model is given in Table 2.

Table 2: An overview of the properties of the models described by Grall et al. (2002); Castanier et al. (2003, 2005) and our new models.

Durations of main-tenance actions Repairs Number of components Identical components

Grall et al. (2002) No No One

-Castanier et al. (2003) Yes Yes One

-Castanier et al. (2005) No No Two, series No

Model 1 No No Two, parallel Yes

(15)

2.2 Maintenance Strategy 9

Here Model 1 and Model 2 are the models that we think are interesting to consider in addition to the original models. Note that besides the properties that are mentioned in Table 2 all models have an aperiodic inspection schedule, a discrete time axis and a fixed failure level. In almost all models failures are detected upon inspection, although in the discrete decision model of Grall et al. (2002) it is assumed that failures are detected immediately. Here a corrective replacement is scheduled at the first inspection time after the moment of failure, instead of at the moment of the next planned inspection. As follows from Table 2 the additional models that we consider have identical components, whereas this does not hold for the model of Castanier et al. (2005). We are adapting this model, which already has quite extensive calculations, so including identical components (with identical thresholds) is a good way to save a lot of computing time. Furthermore we think that it is realistic to consider a system of identical components, think for example of a field with a number of identical pumps. In Table 2 it is also stated that Model 1 and Model 2 have two components that are ‘parallel’. With this we mean that if one component fails the other components keep functioning. In a series system it is assumed that as soon as one component fails the complete system is shut down.

2.2.1 Maintenance Actions

As we mentioned before, both inspections and maintenance actions can only take place at the discrete points in time k. At those moments it is possible, but not mandatory, to monitor the condition or to perform a maintenance operation. Depending on the model that is considered several maintenance actions can be applied. Below these actions are described in somewhat more detail.

Inspection

Determining the current state of the system, and comparing the obtained value with predetermined thresholds.

Corrective replacement

Changing the state of the system to the ‘as good as new’ state after a failure has occurred.

Preventive (condition-based) replacement

Changing the state of the system to the ‘as good as new’ state before a failure occurs. Preventive replacement is initiated when the condition of the component reaches a certain threshold.

Repair

Improving the current state of the system to a better one, not necessarily as good as new. The extent to which the system is improved may be stochastic, and the repair is initiated as soon as the deterioration level exceeds a certain threshold value.

Opportunistic replacement (Two or more components)

Performing a preventive replacement on a component when at least one other component is replaced (either correctively or preventively). In order to be opportunistically replaced the deterioration level of the concerning item should exceed a certain threshold value, while at least one other component should require a replacement.

Opportunistic repair (Two or more components)

(16)

opportunistically repaired the deterioration level of the concerning item should exceed a certain threshold value, while at least one other component should require a repair. (Very similar to opportunistic replacement.)

The reason for including opportunistic replacements and opportunistic repairs is the fact that there can be shared set-up costs involved for repairing or replacing a component. In case several components are maintained at the same time these costs have to be paid only once, which makes it attractive to apply opportunistic maintenance. In practice shared set-up costs can arise from traveling to the right location, scheduling personnel, ordering spare parts or doing paperwork. It is quite plausible that the set-up costs for a replacement of the component are higher than the set-up costs for a repair. Therefore we distinguish between these set-up costs, and we do not allow for an opportunistic repair in case at least one component is replaced but none repaired and vice versa. The resulting costs of the different maintenance operations are given in Table 3. In case of a repair these costs may depend on the improvement in the component’s condition x0.

Table 3: An overview of the costs associated with the different maintenance operations.

Maintenance action Corresponding costs

Inspection (all items) ci

Corrective replacement (per item) cc

Preventive replacement (per item) cp

Repair (per item) cr(x0)

Shared set-up cost repair cs1

Shared set-up cost replacement cs2

When failures are only noticed upon inspection it can occur that the system is not functioning properly for a while before it gets replaced. We therefore introduce a so-called unavailability cost rate cf, incurred per unit of time a component is unavailable due to failure before it is

inspected and replaced. Note that the amount of time a component is unavailable cannot exactly be measured due to the fact that failures are only noticed upon inspection. A possible solution for this is to use an upper bound for the unavailability time. Note that the way these costs are defined here differs somewhat from Grall et al. (2002) and Castanier et al. (2003, 2005), where it is assumed that the inspection costs are included in the repair and replacement costs. Here we do not make this assumption, since it limits the possibilities of extending the model to three or more units.

2.2.2 Multi-Threshold Maintenance Policy

In the articles that we consider the authors propose to use a multi-threshold maintenance policy. The global idea underlying this approach is to determine a set of threshold values for each component on which the scheduling of maintenance activities is based. To that extent we define ξ₀i, ξi₁, . . . , ξ_Ni

i, ξ

i

Ni+1 (with ξ

i

0 ≤ ξ1i ≤ . . . ≤ ξNii+1 ≤ Li) to be the threshold values

of component i, where Ni denotes the fixed number of inspection threshold values for

(17)

whether a replacement or repair should be performed for each item and when to perform the next inspection. For systems consisting of more than one component it is common to in-spect all components at the same time, and for tractability we have to assume such complete system inspections. Below we explain our policy in more detail. First we consider a system which consists of one unit, after which we can extend this policy to the case with two or more components.

One component:

Consider the state of the component (i), by comparing the observed deterioration level with the threshold values.

If ξi_l ≤ Xi

k< ξil+1 for l ∈ {0, Ni− 1}:

Do not perform any maintenance action on component i right now, but schedule the next inspection at Ni− l periods from now.

If ξi_N

i ≤ X

i

k< ξNii+1:

Perform a preventive repair on component i at this moment, which decreases the de-terioration level either stochastically to a value of yi between 0 and ζi, the so-called

restore threshold, or to the fixed amount of yi = ζi depending on the chosen strategy.

If yi ∈ [ξli, ξl+1i ), with l ∈ {0, . . . , Ni − 1}, the next inspection is scheduled at Ni − l

periods from now. If ξi_N

i+1≤ X

i k< Li:

Perform a preventive replacement on component i at this moment, and schedule the next inspection Ni periods later.

If X_ki ≥ Li:

Component i has failed. Perform a corrective replacement on component i at this moment, and schedule the next inspection Ni periods later.

An advantage of this method is that the length of the inspection intervals depends on the de-terioration level of the component. In other words, if a component is relatively new it will not be inspected as often as when it is approaching failure. This can save a lot of time and money. Two or more components:

In case the system consists of two or more components, one can use the same policy as in the case of one component as a first step. In this way it will be determined for each component separately whether the component requires maintenance and when its next inspection should be scheduled. If at least one component needs to be repaired or replaced, the set-up costs need to be paid, and it may be profitable to have other components being repaired or replaced as well. In order to determine whether this is the case we define two additional threshold values, ζ₁i and ζ₂i (with 0 ≤ ζ₁i ≤ ξi

Ni and 0 ≤ ζ

i

2 ≤ ξNii+1). They are used to determine

whether a component should be opportunistically repaired or replaced respectively, given that at least one other component will be repaired or replaced. This can be done by comparing the deterioration level of each other component with the opportunistic threshold values as follows.

If 0 ≤ X_ki < ζ₁i and 0 ≤ X_ki < ζ₂i:

(18)

12 2 SYSTEM DESCRIPTION If ζ₁i ≤ Xi k< ξNii and X j k∈ [ξ j Nj, ξ j

Nj+1) (j repaired) for some j 6= i:

Perform an opportunistic repair on component i at this moment, which decreases the deterioration level to a value of yi, with yi ∈ [ξil, ξl+1i ) and l ∈ {0, . . . , Ni−1}. Reschedule

the next inspection to Ni− l periods from now.

If ζ₂i ≤ Xi

k< ξNii+1 and X

j k∈ [ξ

j

Nj+1, ∞) (j replaced) for some j 6= i:

Perform an opportunistic replacement on component i at this moment, and reschedule the next inspection to Ni periods from now.

Since we assumed that all components need to be inspected at the same time, we need to decide on the next common inspection time. In the previous two steps we have determined the next inspection time for each component separately, which we can now use to find the moment at which the next general inspection will take place according to the following step. Our maintenance policy is developed in such a way that as the deterioration level of a component increases the next inspection is scheduled sooner, so it is important to take into account the components with the highest level of degradation. When considering a system with a small number of components, such as in this Master’s Thesis, it is convenient to use the next inspection time of the components that occurs first, i.e., the minimum of the determined inspection times. At this point in time all components will be inspected. However, for a system consisting of a large number of components this approach may lead to high costs due to the increased number of inspections. In this case a different strategy should be chosen. Example In order to demonstrate our maintenance policy we first consider a system con-sisting of one component. Assume that for this system four maintenance actions are available; inspection, repair, preventive replacement and corrective replacement, all of which take a neg-ligible amount of time. Suppose that we are now at time k and an inspection needs to be performed. Depending on the deterioration level of the component it will be decided whether or not to repair or replace the system and when to perform the next inspection, according to the policy that we explained above. This decision is depicted in Figure 3.

(19)

Figure 3 explains how to compare the observed deterioration level Xk with the threshold

values such that a decision can be made on whether or not to repair or replace the com-ponent and when to schedule the next inspection. It follows that for the lowest values of deterioration, so between ξ0 and ξ1, no further maintenance action is performed and the next

inspection is scheduled N periods later. If the deterioration level is a bit higher, the next in-spection is scheduled somewhat sooner. However, once a deterioration value of ξN or larger is

obtained a maintenance action will be performed before determining the next inspection time. If we assume five threshold values, i.e., three inspection thresholds (N = 3), we can construct a sample path. This path is given in Figure 4. Here we assume that if the component is repaired its deterioration level will be decreased stochastically to a value between zero and its restore threshold ζ.

Figure 4: Sample path of a maintained system with one component, N = 3.

In Figure 4 an example of the wear pattern of a system consisting of one component is given. Just like in Figure 2 the straight lines represent a simplified version of the wear pattern. At time k the component is as good as new, so it has a deterioration level of zero. Since zero is smaller than the threshold value ξ1 the next inspection is scheduled for three periods from

now. At that moment, at k +3, we observe a deterioration level of just below ξ2, which implies

that no repair or replacement is required and that the next inspection is scheduled in two periods. The same holds for k + 5, although now the next inspection is scheduled one period ahead due to the increased deterioration level. We see that the component is repaired at time k + 10, preventively replaced at time k + 6 and correctively replaced at time k + 15. After each replacement the next inspection is scheduled three periods later. In this specific example this also holds for each repair, since the restore threshold ζ is smaller than ξ1. Therefore,

the component’s deterioration level will automatically be between zero and ξ1 after a repair,

(20)

For a system consisting of two or more components it is somewhat more complicated to graphically depict the maintenance policy structure. Nevertheless we would like to give an example path of a two-unit system in order to show the ideas of opportunistic repairs and replacements. Suppose we have a system consisting of two identical components. This means that both components have the same degradation model and can use the same threshold values. An example path of the deterioration patterns of the components is given in Figure 5, where we used the same number of threshold values as in the previous example (N = 3).

(21)

15

3 Mathematical Model

In this section we construct the mathematical expressions for the models that we consider. We consider multiple components, as the case with one component can be considered as a special case. We start off by explaining that we can view our models as semi-regeneration cycles, which means that we can construct a stationary law of the maintained system state. Afterwards we construct the performance criteria, for which expressions can be derived by using this law. We mathematically formulate Model 1, while the formulation of Model 2 can be found in the Appendix.

3.1 Semi-Regeneration Cycles and Stationary Law

In the previous section it is explained that after a replacement of the complete system, either preventively or correctively, the system state is as good as new and its deterioration does not depend on the past. This implies that we can view our deterioration process (Xk)k∈N

as a regeneration process, of which a proof is given in Grall et al. (2002). The regeneration points are given by the points in time the system has just been replaced and, hence, is as good as new. This means that in order to evaluate our system we can consider the long run behavior of the process in a single renewal cycle, i.e., the period between two successive system replacements. Although this simplifies calculations significantly, especially for multiple components the model is still hard to solve. Grall et al. (2002) and Castanier et al. (2003, 2005) explain that we can go one step further and observe the period between two successive system inspections as semi-regeneration cycles, since the decision to be made at the moment of an inspection solely depends on the system state at that point in time. A drawback of considering semi-regeneration cycles is that it imposes some restrictions on the model, e.g., when including durations of maintenance activities. In Figure 6 the example path of Figure 4 is given, along with the regeneration cycles R and semi-regeneration cycles S.

(22)

16 3 MATHEMATICAL MODEL

The calculations for, for example, the long run average costs of maintaining a system can be based on the average costs in one inspection cycle in steady-state. The existence of semi-regeneration cycles means that we can construct a probabilistic law π(x1, x2, . . . , xn)

to describe the probability of being in a certain state at the moment of an inspection. The basic idea of the used approach is that the probability of being in state x at this inspection is equal to the probability of being in state y at the previous inspection multiplied by the probability of moving from state y to state x during the semi-regeneration cycle, aggregated over all possible values of y. Of course the exact stationary laws of the various models that we consider differ somewhat from each other, which is why we show the exact calculations for each model separately. The logic is to first consider the possible steps during an inspection cycle, i.e., performing a repair, a replacement, or no maintenance action at all, and then to use these possibilities to calculate the stationary law.

3.2 Performance Criterion

As a performance criterion for our models we use the long run average operating cost per unit of time. This is also done in Grall et al. (2002) and Castanier et al. (2005), while in Castanier et al. (2003) in addition the long run average availability is considered. The operating costs consists of all costs involved in performing inspections, repairs, preventive replacements, corrective replacements, the set-up costs for repairs and replacements, and the system unavailability due to a failure summed over all components. Let C(t) denote the cumulative operating cost up to time t, and let CI(t) denote the inspection costs of all

components, C_Ri (t), C_Pi(t) and C_Ci(t) denote the costs resulting from repairing, preventively replacing and correctively replacing component i respectively (not including the set-up costs), Ci

F(t) denote the costs for unavailability of component i due to a failure, and CS1(t) and

CS2(t) respectively denote the set-up costs for performing repairs or replacements. All costs

are cumulated up to time t. It follows that the cumulative operating costs C(t) are constructed as follows. C(t) = CI(t) + n X i=1 C_Ri(t) + C_Pi(t) + C_Ci(t) + C_Fi(t) +CS1(t) + CS2(t) (1)

Note that depending on the chosen maintenance strategy some components of this cost equa-tion will be left out, for example if one does not consider repairs.

In order to obtain the long run average operating costs we can make use of Equation (1) as follows. Let C∞ denote the long run average operating costs. Then

C∞ = lim

t→∞

E[C(t)]

t .

However, in the previous subsection we explained that we can view the periods between two successive inspections as semi-regeneration cycles. This means that, at steady state, we can write C∞ as follows, with S the length of an inspection cycle.

C∞ =

Eπ[C(S)]

Eπ[S]

(23)

3.3 Calculations Per Model 17

3.3 Calculations Per Model

In this part we will show the calculations for Model 1, while the calculations for Model 2 are given in the Appendix. Since the models of Grall et al. (2002) and Castanier et al. (2003, 2005) are already described in the respective articles, we will not repeat those equations. Afterwards, we explain some difficulties that arise from including durations of maintenance activities in Model 1 and from extending this model to include three components.

3.3.1 Model 1

Model 1 is similar to the model in Castanier et al. (2005), but instead of a two-unit series system we consider two units in parallel. This means that if one of the components fails the other can still function, whereas in the mentioned article the whole system fails as soon as one of the components fails. Furthermore we consider two identical components, which implies that we can use identical distribution functions and threshold values for both components. I.e., N1= N2= N , L1 = L2 = L and ξi1= ξi2= ξi for i = 0, 1, . . . , N + 1. We do not consider

repairs, so we can set the repair threshold ξN equal to ξN +1to prevent repairs from occurring,

along with ζ1= ξN +1.

Cost Components We can write the components of the cost function (2) as follows. First, some components of the cost function are equal to zero, since we do not consider repairs.

Eπ[CRi (S)] = 0, i = 1, 2 (No repairs)

Eπ[CS1(S)] = 0 (No repairs)

Second, we know that in each inspection cycle the inspection costs will have to be paid once. Hence,

Eπ[CI(S)] = ci.

The preventive replacement costs cpfor component i arise either when component i requires a

preventive replacement, independent of the condition of component j, or when component j is replaced (preventively or correctively) and component i can be opportunistically replaced. We can express this in the following way. Note that the set-up costs will be calculated separately.

Eπ[CPi(S)] = cp·      Z L ξN +1 Z ∞ 0 π(x1, x2)dxjdxi | {z } i preventively replaced + Z ξN +1 ζ2 Z ∞ ξN +1 π(x1, x2)dxjdxi | {z } i opportunistically replaced      , i = 1, 2

The corrective replacement costs ccfor component i are incurred when component i has failed,

independent of the condition of component j. This implies that Eπ[CCi (S)] = cc· Z ∞ L Z ∞ 0 π(x1, x2)dxjdxi, i = 1, 2.

(24)

costs to be paid is given by

Eπ[CS2(S)] = cs2· P(One or more replacements)

= cs2·      Z ∞ ξN +1 Z ζ2 0 π(x1, x2)dx2dx1 | {z } Component 1 replaced + Z ζ2 0 Z ∞ ξN +1 π(x1, x2)dx2dx1 | {z } Component 2 replaced + Z ∞ ξN +1 Z ∞ ζ2 π(x1, x2)dx2dx1+ Z ξN +1 ζ2 Z ∞ ξN +1 π(x1, x2)dx2dx1 | {z } Both replaced      = cs2· Z ∞ ξN +1 Z ∞ 0 π(x1, x2)dx2dx1+ Z ξN +1 0 Z ∞ ξN +1 π(x1, x2)dx2dx1 !

The average length of a semi-regeneration cycle can be calculated as follows. Here we dis-tinguish three cases; no replacements, one replacement (either on component 1 or 2) and a complete system renewal.

Eπ[S] = N X k=1 k · P(S = k) = N X k=1 k ·      Z ξN −k+1 ξN −k Z ξN −k+1 0 π(x1, x2)dx2dx1+ Z ξN −k 0 Z ξN −k+1 ξN −k π(x1, x2)dx2dx1 | {z } No replacements + Z ∞ ξN +1 Z min{ζ2,ξN −k+1} min{ζ2,ξN −k} π(x1, x2)dx2dx1+ Z min{ζ2,ξN −k+1} min{ζ2,ξN −k} Z ∞ ξN +1 π(x1, x2)dx2dx1 | {z } One replacement + I{k=N }· Z ∞ ξN +1 Z ∞ ζ2 π(x1, x2)dx2dx1+ Z ξN +1 ζ2 Z ∞ ξN +1 π(x1, x2)dx2dx1 ! | {z } Two replacements      

Just like in Castanier et al. (2005) it is impossible to trace the exact failure time of a com-ponent due to the discrete time horizon we consider. Therefore, the time spent in the failure zone D_Fi (t) during time t for component i, i=1,2, needs to be estimated. Suppose we consider the time span [0, S], and assume that a component fails at time tf, somewhere between times

k and k + 1. Since the levels of deterioration at k and k + 1 are known, a possible way to esti-mate the failure time of a component is to assume a linear increase in deterioration between these two moments. This leads to an estimated moment of failure of L−Xk

Xk+1−Xk time units later

than k, so the failure time of the component can be estimated by S −k + L−Xk

Xk+1−Xk

(25)

However, in Castanier et al. (2003, 2005) the failure time of a component is estimated by an upper bound ¯Di_F(t). Suppose we consider the same time span [0, S] and moment of failure tf as before. The upper bound is then given by the longest possible time spent in the failed

state, which is equal to S − k, as shown in Figure 7. Here tf denotes the moment of failure.

Figure 7: Exact failure time Di_F(S) of component i and its upper bound ¯D_Fi (S). It follows that the difference between the estimated failure time for our intuitive method and the upper bound as used in Castanier et al. (2003, 2005) is equal to at most one time unit. We decide to use the latter method, although the corresponding costs will be somewhat higher than with our own method, since this enables us to compare our results with those of Castanier et al. (2005). Besides, it is easy to change this assumption in a later stage.

We can now express the unavailability cost of each component as a function of this upper bound as follows. Eπ[CFi (S)] = cf· Eπ[DFi (S)] ≈ cf· Eπ[ ¯DFi (S)] = cf· N X k=1 k · P( ¯Di_F(S) = k)

In order to calculate the function P( ¯Di_F(S) = k) we first define Hi(k|y; l) as the probability

that the unavailability time of component i is equal to k time units, given that after the previous inspection and possible maintenance actions the component had a deterioration level of y and the next inspection was scheduled l time units later, i.e., ξN −l ≤ y < ξN −l+1.

Then Hi(k|y; l) is given by

Hi(k|y; l) =        0 l < k R∞ L−yf (t)dt l = k RL−y 0 f (l−k)_(u)R∞ L−y−uf (s)ds du l > k , y ∈ [ξN −l, ξN −l+1).

We can use this function in order to further calculate Eπ[C_Fi(S)], for which we distinguish four

cases; no replacement, only component i being replaced, only component j being replaced (j 6= i) and a complete system replacement.

Eπ[CFi(S)] ≈ cf · N

X

k=1

(26)

20 3 MATHEMATICAL MODEL = cf · N X k=1 k ·        N X li=1 N X lj=1 Z ξ_{N −li+1} ξ_{N −li} Z ξ_{N −lj +1} ξ_{N −lj}

π(y1, y2)dyjHi(k|yi; min{li, lj})dyi

| {z } No replacement + N X lj=1 Z ∞ ξN +1 Z min{ζ2,ξ_{N −lj +1}} min{ζ2,ξ_{N −lj}} π(y1, y2)dyjdyiHi(k|0; lj) | {z } i replaced + N X li=1 Z min{ζ2,ξ_{N −li+1}} min{ζ2,ξ_{N −li}} Z ∞ ξN +1

π(y1, y2)dyjHi(k|yi; li)dyi

| {z } j (j 6= i) replaced + Z ∞ ξN +1 Z ∞ ζ2 π(y1, y2)dyjdyi+ Z ξN +1 ζ2 Z ∞ ξN +1 π(y1, y2)dyjdyi ! · H_i(k|0; N ) | {z } Both replaced      

Probability Law In order to actually calculate the cost function (2) the probability law π(x1, x2) is needed. As we explained before, π(x1, x2) denotes the probability density function

of being in state (x1, x2) at the start of an inspection. It can be calculated as the

probabil-ity of being in state (y1, y2) multiplied by the probability of moving from state (y1, y2) to

state (x1, x2) during an inspection cycle. Here we distinguish three cases; no replacement,

one component being replaced and a complete system renewal initiated after the previous inspection.

No replacement

If y1 ∈ [ξk, ξk+1) and y2 ∈ [ξl, ξl+1) with k, l ∈ {0, 1, . . . , N − 1} no replacement will

be performed. The next inspection will be scheduled min{N − k, N − l} = N + min{−k, −l} = N −max{k, l} periods later. The increases in deterioration levels from y1

to x1and from y2to x2have a density of f(N −max{k,l})(x1−y1) and f(N −max{k,l})(x2−y2),

respectively.

One replacement

If yi ∈ [ξN +1, ∞) and yj ∈ [ξk, ξk+1) with yj < ζ2 (j 6= i) only component i will be

replaced. The next inspection is scheduled N − k periods later, while component i gets a deterioration level of 0 due to the replacement. During the N − k periods before the next inspection the increases in deterioration levels from 0 to xi and from yj to xj have

a density of f(N −k)(xi) and f(N −k)(xj − yj) respectively. This holds for i = 1, 2.

Two replacements

If y1 ∈ [ξN +1, ∞) and y2 ∈ [ζ2, ∞) or if y1 ∈ [ζ2, ξN +1) and y2 ∈ [ξN +1, ∞) both

components are replaced. The next inspection is scheduled N time units later, which implies that the increases in deterioration levels from 0 to x1 and from 0 to x2 have a

(27)

If we use this information we can construct the probability law π(x1, x2) as follows.

π(x1, x2) = N −1 X k=0 N −1 X l=0 Z ξk+1 ξk Z ξl+1 ξl

π(y1, y2)f(N −max{k,l})(x1− y1)f(N −max{k,l})(x2− y2)dy2dy1

| {z } No replacement +X i6=j N −1 X k=0 Z min{ζ2,ξk+1} min{ζ2,ξk} Z ∞ ξN +1 π(y1, y2)dyi ! f(N −k)(xj− yj)dyj ! f(N −k)(xi) | {z } One replacement + Z ∞ ξN +1 Z ∞ ζ2 π(y1, y2)dy2dy1+ Z ξN +1 ζ2 Z ∞ ξN +1 π(y1, y2)dy2dy1 ! f(N )(x1)f(N )(x2) | {z } Two replacements

3.3.2 Including Durations of Maintenance Activities

(28)

3.3.3 Including Three Components

(29)

23

4 Performance Assessments of the Maintenance Policies

In the previous section we derived mathematical expressions for Model 1, which we can use in this section to asses the performances of this model, after programming it in R (R Core Team (2012)). In order to do so we first provide some underlying theory on integral equations, successive approximation and numeric integration. Afterwards we rewrite our probability law as a Fredholm integral equation, which is required in order to be able to solve it. Then we assess the performances of the model.

4.1 Underlying Theory

In Castanier et al. (2005) it is stated that the probability law, which is a bi-dimensional one-sided integral equation of the second kind, can be solved by using an iterative method based on the method of successive approximations, due to the regularity of the densities f_i(n). This is confirmed in Castanier et al. (2003), in which it is explained that the probability law can be numerically evaluated by using an iterative algorithm adapted from classical integration schemes for Volterra equations. We will first give a short summary on different types of integral equations, such as Volterra and Fredholm equations. Afterwards we explain the method of successive approximation. Since we will solve the integral equations numerically, we also explain the method of numeric integration that we use.

4.1.1 Integral Equations

In Castanier et al. (2003) it is stated that in order to solve the probability law corresponding to their model an iterative algorithm adapted from classical integration schemes for Volterra equations can be used. In this section we give a short introduction on different types of Fredholm and Volterra integral equations.

In Polyanin and Manzhirov (1998) and Press, Teukolsky, Vetterling, and Flannery (1992) an overview of different types of integral equations is given, which provides the basis for this overview. Suppose we consider the following integral equation.

βy(x) = F (x) + Z b

a

K(x, t)y(t)dt (3)

Here x is the independent variable, y(x) the unknown function, K(x, t) is called the kernel, and F (x) is some known function of x. Equation (3) is the general expression for a Fredholm integral equation. In case the kernel K(x, t) = 0 for t > x this equation is called a Volterra integral equation, and it can be written as follows.

βy(x) = F (x) + Z x

a

K(x, t)y(t)dt (4)

(30)

24 4 PERFORMANCE ASSESSMENTS OF THE MAINTENANCE POLICIES

4.1.2 Successive Approximation

We are forced to use a numeric method for obtaining an expression for the probability law due to the its complexity. In Castanier et al. (2005) the method of successive approximation is used for this purpose. Here we explain this method. Suppose we consider a system consisting of one unit (e.g. the model developed in Castanier et al. (2003). The probability law that we need for that model, which is currently specified as a homogeneous Fredholm integral equation of the second kind, can be written as a nonhomogeneous Fredholm integral equation of the second kind, i.e.,

π(x) = F (x) + Z b

a

K(x, t)π(t)dt. (5)

This reformulation is required for applying the method of successive approximation. Note that in Section 4.2 we rewrite the probability law of Model 1 in this way, while the reformulation of Model 2 is given in the Appendix. We demonstrate the logic of the method of successive approximation by rewriting Equation (5) as follows.

(31)

4.1 Underlying Theory 25

From this it follows that we can calculate the probability law as π(x) = lim

n→∞πn(x)

Provided that the series converges, we can now approximate the probability law by obtaining an expression for πn(x) for some sufficiently large n. This procedure is called successive

approximation. In Castanier et al. (2005) it is explained that this method can be applied due to the regularity of the density functions fi. Note that in our policies as described in Model 1

and Model 2 we consider two-dimensional probability laws. In these cases the method of successive approximation can also be applied.

4.1.3 Numeric Integration

Due to the high level of complexity of the equations that we need to solve we are forced to use numeric integration. In Press et al. (1992) it is explained that the extended midpoint rule is a classical formula used to solve integrals numerically, which is still rated well. Moreover, computing time issues prevent us from evaluating the function to be integrated at a lot of different points, so we are forced to use a simple rule such as this one.

Extended midpoint rule Below we explain the extended midpoint rule, first for the case with a finite upper bound, and afterwards for the case with an infinite upper bound.

Finite Upper Bound Suppose that we would like to integrate the function f between a lower limit a and an upper limit b. The extended midpoint rule is an open formula, since it does not require evaluation of the function value at the bounds of integration. Instead, the function is divided in a number M (M ≥ 1) of subintervals for which the function value is evaluated at the middle of the interval and used for approximating the integral in that subinterval. This procedure is demonstrated in Figure 8, where the extended midpoint rule is applied to the function f (x) = x2_{− 2x. Here we integrate from zero to four, and divide}

this interval into four subintervals of equal length.

(32)

In case M = 1 we have the midpoint rule. In Antia (2002) it is explained that this rule states that Z b a f (x)dx = (b − a) · f b − a 2 + O (b − a) 3_f00 24

where O(·) denotes the error term, which depends on b − a and on the second derivative of the function f evaluated somewhere in the interval between a and b. As its name suggests, the extended midpoint rule extends this rule by allowing for more subintervals (M ≥ 1) to approximate the integral. In Press et al. (1992) it is explained that the extended midpoint rule is given by Z b a f (x)dx = b − a M · M X i=1 f a + (i −1 2) b − a M + O 1 M2

Infinite Upper Bound Instead of integrating between a lower limit a and a finite upper limit b, suppose now that we would like to integrate the function f between a lower limit a and infinity. From the definition of the extended midpoint rule above it immediately follows that this rule cannot be applied in this case, since it would require an infinite number of steps. Changing the variables such that we would integrate on a finite interval solves this problem, due to the open form of the formula. To illustrate this, we rewrite our integral equation as follows. Z ∞ a f (x)dx = Z 1 0 f a + x 1 − x 1 (1 − x)2dx (7)

To show that this is true we can simplify the right-hand side of the equation by using substitu-tion (Sydsaeter and Hammond (2008)). Let u = a +_1−xx . Then it follows that du = _(1−x)1 2dx.

The bounds of integration 0 and 1 will be replaced by a + ₁₋₀0 = a and lim

t→1 a + t

1−t = ∞.

Hence, we can rewrite the right-hand side of Equation (7) as Z 1 0 f a + x 1 − x 1 (1 − x)2dx = Z ∞ a f (u)du,

from which it follows that Equation (7) is correct. After rewriting the equation we can simply apply the extended midpoint rule as explained above.

An alternative solution for dealing with infinite upper bounds is to use some finite upper bound instead. In our maintenance policy we assume that a component i fails as soon as its deterioration level exceeds the failure level Li. It will then be replaced at the next planned

inspection, after which the component is as good as new. Taking this into account, it makes sense to assume that the deterioration level of a component will never exceed a certain level, for example 2.5 times the failure level Li.

Note that the midpoint rule converges rather slowly, about 200 function evaluations are re-quired to get an accuracy of seven decimal figures (Antia (2002)).

(33)

4.2 Rewriting the Probability Law 27

to the following. Suppose that we immediately apply the logic of the extended midpoint rule to the two-dimensional case. This means that we divide the area to be integrated in M2 different (sub)squares. The integral in each square is then approximated by evaluating the function value at the middle of the corresponding square and multiplying it with the area of the square. We can depict this situation as in Figure 9.

Figure 9: Extended midpoint rule applied to f (x1, x2), integrated from a to b and from c to

d, with M = 4.

We will use this method for numerically integrating two-dimensional equations, along with M = 25. Due to computing time issues we cannot select a larger value for M .

4.2 Rewriting the Probability Law

When we described the mathematical model of the maintenance policy that we consider in Model 1, we also gave an expression for the corresponding probability law. Currently the probability law is formulated as a homogeneous Fredholm equation of the second kind, which means that F (x) equals zero in Equation (5). This implies that if we now apply the method of successive approximation all estimates will equal zero. In order to avoid this problem we first rewrite our probability law to a nonhomogeneous Fredholm equation of the second kind. Note that the reformulation of Model 2 is given in the Appendix.

4.2.1 Model 1

In Section 2 we formulate the probability law of Model 1 as follows. π(x1, x2) = N −1 X k=0 N −1 X l=0 Z ξ_k+1 ξk Z ξ_l+1 ξl

(34)

28 4 PERFORMANCE ASSESSMENTS OF THE MAINTENANCE POLICIES +X i6=j N −1 X k=0 Z min{ζ2,ξk+1} min{ζ2,ξk} Z ∞ ξN +1 π(y1, y2)dyi ! f(N −k)(xj− yj)dyj ! f(N −k)(xi) + Z ∞ ξN +1 Z ∞ ζ2 π(y1, y2)dy2dy1+ Z ξN +1 ζ2 Z ∞ ξN +1 π(y1, y2)dy2dy1 ! f(N )(x1)f(N )(x2)

Note that we can rewrite this function as π(x1, x2) = Z ∞ 0 Z ∞ 0 π(y1, y2) · K(x1, x2, y1, y2)dy2dy1

by specifying the kernel K(x1, x2, y1, y2) as follows.

K(x1, x2, y1, y2) =                                    f(N )(x1)f(N )(x2) If y1∈ [ξN, ∞) ∩ y2∈ [ζ2, ∞) or y1∈ [ζ2, ξN) ∩ y2∈ [ξN, ∞)

f(N −i)(x1− y1)f(N −i)(x2) If y1∈ [min{ξi, ζ2}, min{ξi+1, ζ2})∩

y2∈ [ξN, ∞), i = 0, 1, . . . , N − 1

f(N −i)(x1)f(N −i)(x2− y2) If y2∈ [min{ξi, ζ2}, min{ξi+1, ζ2})∩

y1∈ [ξN, ∞), i = 0, 1, . . . , N − 1

f(N −i)(x1− y1)f(N −i)(x2− y2) If y1∈ [0, ξi) ∩ y2∈ [ξi, ξi+1) or

y1∈ [ξi, ξi+1) ∩ y2 ∈ [0, ξi+1),

i = 0, 1, . . . , N − 1

Now the probability law is written as a homogeneous Fredholm integral equation of the second kind. We notice that the kernel currently contains a part that does not depend on y1 or y2.

Therefore, it is possible to simplify the integral equation somewhat, as we show below. In our calculations we make use of the fact that R₀∞R₀∞π(x1, x2)dx2dx1= 1, since π is a probability

(35)

4.3 Performance Assessments 29 + 1 − Z ζ2 0 Z ∞ 0 π(y1, y2)dy2dy1+ Z ξN ζ2 Z ξN 0 π(y1, y2)dy2dy1 + Z ∞ ξN Z ζ2 0 π(y1, y2)dy2dy1 · f(N )(x1)f(N )(x2) = F (x1, x2) + Z Z S π(y1, y2) ¯K(x1, x2, y1, y2)dy2dy2

which is the specification of a nonhomogeneous Fredholm integral equation, with F (x1, x2) = f(N )(x1)f(N )(x2)

S = {(y1, y2) : (y1∈ [0, ξN) ∩ y2∈ [0, ξN)) ∪ (y1 ∈ [0, ζ2) ∩ y2 ∈ [ξN, ∞))

∪ (y1 ∈ [ξN, ∞) ∩ y2 ∈ [0, ζ2))}

¯

K(x1, x2, y1, y2) = K(x1, x2, y1, y2) − f(N )(x1)f(N )(x2) (8)

Now the method of successive approximation can be applied to the probability law in order to obtain an expression for it.

4.3 Performance Assessments

In this section we assess the performances of Model 1. In order to calculate the expected long run costs per maintenance period we need that the probability law converges in applying the method of successive approximation. This will be investigated in this section. Next we describe a case that we will consider, by specifying the required distribution functions and the different costs of the maintenance actions, along with the used number of thresholds. For this specific case we determine the (near) optimal solution, and explore the influence of the threshold values. We consider what happens if we increase the number of thresholds, and we compare the solution to some classical maintenance policies. Afterwards we consider what happens if we decrease the unavailability cost rate and if we increase the set-up costs for replacements.

4.3.1 Convergence of Probability Law

As we explained in Section 4.1.2 we apply the method of successive approximation to obtain our probability law π(x1, x2). Note that this probability law depends on the values of the

thresholds. Suppose we consider a policy with N = 2, the inspection threshold ξ1 = 1.4,

the preventive replacement threshold ξ2 = 1.6, and the opportunistic replacement threshold

ζ2 = 1.1. Here ξ2 denotes the level at which a preventive replacement will be performed.

Furthermore, we assume that the failure level is L = 2 and that the increase in deterioration level per time unit follows an exponential distribution with parameter α = 3. For this specific case we determine πn(x1, x2) for n = 0, 1, . . . , 40 to see how fast it converges. For example

(36)

Figure 10: πn(0.5, 3.3) for n = 0, 1, . . . , 40.

From Figure 10 it appears that the probability law converges rather quickly. Already after seven iterations the probability law does not seem to change anymore. For other values of the thresholds we either find the exact same result, or we see that the probability law diverges rapidly. The latter case can be caused by the fact that for some values of x1, x2, y1 and y2

the function ¯K(·) (in Equation (8)) is negative. Although this should not cause any trouble, the approximation of its integral can yield negative outcomes due to the fact that we only use M = 25 in the numeric integration. This, in turn, can cause some values of πnto be negative.

In practice it is of course not possible for a distribution function to have negative values, so we decided to ignore the corresponding solution if this is the case.

Based on the previous, we decide to use 15 iterations for our calculations of the average long run costs per maintenance period, and to verify whether the probability law does not contain negative values before using it. For the case specified above, i.e., L = 2, α = 3, ξ1 = 1.4, ξ2=

1.6 and ζ2 = 1.1, this leads to the probability law shown in Figures 11 and 12, where Figure

12 uses colors to represent the value of π(x1, x2). Note that the probability law is symmetric

with respect to the line x1= x2. This can be explained by the fact that we are dealing with

identical components. Furthermore, we assumed that the deterioration levels will not exceed the value of 2.5 times L, and we observe that for x1 and x2 approaching five π(x1, x2) indeed

has a value which is very close to zero. From Figure 12 we observe that the maximum value of π is obtained for x1 = x2 = 1.5, which is located between the threshold values ξ1 and

ξ2. In case either x1 or x2 is approximately equal to 1.5 we also observe high values of π.

This can be explained by the following. If, after performing maintenance actions, at least one of the deterioration levels exceeds ξ1 = 1.4 the next inspection is scheduled one period

later, otherwise two periods. This explains the larger probability mass for x1 or x2 exceeding

ξ1. However, a component will be preventively replaced when an inspection reveals that its

deterioration level exceeds the value of ξ2 = 1.6, and correctively replaced when it exceeds

(37)

4.3 Performance Assessments 31

Figure 11: Probability law π(x1, x2) for ξ1= 1.4, ξ2= 1.6 and ζ2 = 1.1.

(38)

4.3.2 Scenario under Investigation

For determining a case for which we can assess the performances we consider Castanier et al. (2003, 2005). For the performance assessments of their policies these authors choose to use (among others) the costs shown in Table 4.

Table 4: Costs of Example Scenario

Maintenance action Corresponding costs

Inspection ci= 1

Preventive replacement cp= 40

Corrective replacement cc= 100

Unavailability (per time unit) cf = 10

Set-up for replacements cs2 = 20

However, if we compare the unavailability cost rate with the costs corresponding to a pre-ventive and a corrective replacement we observe that cf is rather small. This may indicate

that it is cheaper to let the system be in the failed state than to maintain the components. In practice this is not a realistic case, so we decide to use cf = 1000 instead of cf = 10. We

will not adapt the other costs.

In Castanier et al. (2003, 2005) it is explained that the exponential distribution is a common choice for a distribution function describing the deterioration of a component during one time unit. Hence, f is assumed to be exponential. This implies that the sum of n exponentially distributed variables, or the deterioration gained in n consecutive time units, follows a Gamma distribution. The parameter of the exponential distribution is set to α = 3. Moreover, the failure level is assumed to be L = 2, which means that a component fails as soon as its deterioration level exceeds the value of two. If we consider a component that is as good as new at time zero, and that does not undergo any maintenance actions, its failure probability will evolve as in Figure 13.

(39)

From Figure 13 it follows that N should not be chosen too large, especially when considering relatively high unavailability costs. Note that the average time to a failure is equal to six time units. We decide to set the number of inspection thresholds to N = 2, so the next inspection is planned either one or two periods later. We now need to optimize the inspection thresholds ξ1 and ξ2 as well as the opportunistic replacement threshold ζ2. In R it takes on average

161.89 seconds to calculate the costs corresponding to one such scenario on a laptop with a 1.83GHz dual core processor and 1.99 GB of RAM. In order to give an idea of the behavior of the cost function for different realizations of the threshold values we create Figure 14. In this figure the long run average costs per time unit for the policy are given for different values of ξ2 on the x-axis, of ζ2 on the y-axis, and of ξ1 in different colors.

Figure 14: C∞ for different values of ξ1, ξ2 and ζ2.

We calculated the value of C∞ at ξ1, ξ2, ζ2 ∈ {0, 0.2, 0.4, . . . , L(= 2)}, with 0 ≤ ξ1 ≤ ξ2 and

0 ≤ ζ2 ≤ ξ2, in order to construct Figure 14. Our next observations will be based on this

dataset as well. Note that in Figure 14 there is no value of C∞ for ξ1 > ξ2 and for ζ2 > ξ2,

which means that the upper-right half of the figure is empty. We also observe that for some values of ξ2 and ζ2 close to two the solution is omitted, since the probability law contains

negative values. For large values of ξ2 it holds that ξ1 has a large influence on the value of

C∞, while ζ2 influences the solution to a smaller extent. Overall, the behavior of C∞ seems

(40)

we create plots in which we show the minimal value of C∞ given one of the threshold values.

The result is shown in Figure 15. We also make plots in which we fix one of the threshold values and investigate for which realizations of the other threshold values the minimal costs are obtained. The resulting plots are given in Figure 16.

Figure 15: The minimal value of C∞for

dif-ferent values of ξ1, ξ2 and ζ2.

(41)

With the current set of values for the thresholds at which the performance criterion is evalu-ated we can get a good view on the behavior of the expected long run costs per maintenance period, since from Figure 15 it follows that each graph seems to be rather smooth. We no-tice that the minimum of C∞ is strictly increasing for fixed values of ξ1, the threshold value

that indicates whether the next inspection should be scheduled one time unit later instead of two. Hence, it is profitable to perform inspections quite often, which can be explained by the high unavailability cost rate and the relatively low costs corresponding to an inspection. The threshold value ξ2 indicates whether or not to preventively replace a component. Apparently,

fixing this threshold to a rather high value increases the maintenance costs, due to the increas-ing failure probability. On the other hand, if ξ2 is fixed to a value close to zero maintenance

costs will increase as well, since the preventive replacements will be performed too early in the process. The threshold value that corresponds to opportunistic replacements is given by ζ2.

We notice that C∞ increases for fixed values of ζ2, with ζ2 ≥ 0.6, while it fluctuates around

38 for smaller values of ζ2. This can be caused by the fact that it must hold that ξ2 ≥ ζ2.

For large values of ζ2 less preventive replacements will be performed, which decreases the

availability of the system, and hence increases the costs. For smaller values of ζ2 this problem

does not occur. In Figure 16 we show the values of the thresholds corresponding to the min-imal values of C∞that are shown in Figure 15. If we consider a fixed ξ1 we observe that the

difference between the optimal value of ξ2 and ξ1 is decreasing as ξ1 increases. For smaller

values of ξ1the system is inspected more often than for larger values of ξ1, so the deterioration

level will increase less between two consecutive inspections for small values of ξ1than for large

values of ξ1. Hence, the threshold ξ2 can be increased somewhat when increasing ξ1, until

for ξ1 ≥ 1.2 both thresholds are equal to each other. This is due to the high unavailability

cost rate and the fact that we need ξ2 ≥ ξ1. Furthermore, for fixed values of ξ1 we observe

that the optimal value of ζ2 is first increasing somewhat, after which it decreases to a value

of zero for ξ1≥ 1.8. When increasing ξ1 (and ξ2) less preventive replacements are performed,

so the level at which an opportunistic replacement will be performed is decreasing in order to prevent failures from occurring, while it should be taken into account that ζ2 cannot exceed

ξ2. For fixed values of ξ2 we need that both ξ1≤ ξ2 and ζ2 ≤ ξ2. We observe that for ξ2 ≤ 1.2

the system is inspected quite often, while ζ2 is increasing somewhat. For ξ2 exceeding 1.2

preventive replacements will be performed less regularly, so it becomes more profitable to perform opportunistic replacements. However, the system will be inspected less often, since it will take some time for the deterioration level to reach the high threshold level ξ2. If we fix

the value of ζ2 we observe that when opportunistic replacements are performed quite often,

i.e., when ζ2 is close to zero, the optimal value of ξ2 is increasing somewhat. For ζ2 ≥ 1.2 it

holds that ξ2 = ζ2, and hence no opportunistic replacements will be performed. This can be

explained by the fact that ξ2 should not be chosen too large due to the high unavailability

cost rate, while we need that ξ2 ≥ ζ2. For smaller values of ζ2 we do have opportunistic

re-placements, since the optimal value of ξ2 is exceeding ζ2. Note that also ξ1 cannot exceed ξ2.

We could say that as the opportunistic threshold increases, less inspections will be performed. Results indicate that the minimal average long run cost is located somewhere around 54.344 per period, for ξ1 = 0, ξ2 = 0.6 and ζ2= 0.2. We refer to this solution as the optimal solution,

although we only considered ξ1, ξ2, ζ2 ∈ {0, 0.2, 0.4, . . . , L(= 2)}. It could be that a better

(42)

replaced when its deterioration level exceeds ξ2 = 0.6, and a component will be

opportunis-tically replaced when its deterioration level exceeds ζ2 = 0.2, while the other component has

to be replaced either preventively or correctively. Since ξ1 = 0 the system will be inspected

during each period, except when both component are replaced at the same time. Then the next inspection will be scheduled two periods later. We will now investigate whether we can improve this solution by increasing the number of inspection thresholds from N = 2 to N = 3, and how this result compares to results obtained by classical maintenance policies. Also the cases in which the unavailability costs rate and the set-up cost for a replacement are adapted will be investigated in this section.

4.3.3 Adding a Threshold Value

Let ξ1 = 0, ξ2 = 0.6 and ζ2 = 0.2. For this set of threshold values the expected long run costs

of the maintenance policy are equal to 54.344 per period. We assume that these threshold values are given, and that we add one extra threshold ξnew. There are three possible cases,

either ξnew = 0, 0 < ξnew ≤ 0.6 or 0.6 < ξnew≤ L. The resulting costs for all cases are given

in Figure 17.

Figure 17: C∞ for different values of ξnew (with ξ1= ξ2 = 1.4 and ζ2 = 1).

Note that in Figure 17 the ‘kinks’ of C∞ are located at ξnew = ζ2 = 0.2 and ξnew= ξ2 = 0.6.

It follows that the lowest costs are obtained for ξnew = 0.2, with a long run average cost

of about 35.417 per period. The corresponding set of threshold values is ξ1 = 0, ξ2 = 0.2,

ξ3 = 0.6 and ζ2 = 0.2. This means that the next inspection is scheduled either one, two or

A Condition-Based Maintenance Policy for a Continuously Deteriorating Multi-Unit System with Aperiodic Inspections

A Condition-Based Maintenance Policy for a Continuously

Deteriorating Multi-Unit System with Aperiodic Inspections

A Condition-Based Maintenance Policy for a Continuously

Deteriorating Multi-Unit System with Aperiodic Inspections

Contents

Nomenclature

1

Introduction

2

System Description

3

Mathematical Model

4

Performance Assessments of the Maintenance Policies