# A Generalized Chart-Based Decision-Making Tool for Optimal Preventive Maintenance Time under Perfect Renewal Process Modeling

Authors: Pablo Viveros Gunckel; Fredy Kristjanpoller; René Tapia Peñaloza
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

The most commonly used probabilistic model in reliability studies is the Perfect Renewal Process (PRP), which is characterized by the condition or type of maintenance represented: once the maintenance activities are executed, the equipment is restored to its original condition, leaving it “as good as new.” It is widely used since it represents an optimistic state when an item is replaced, assuming a perfect operational condition of the item after the maintenance. Some models have been developed for determining optimum preventive maintenance (PM) based on different criteria, and almost all aimed at PRP reliability modeling. The contribution of this paper is to analyze a model for determining the optimal preventive maintenance policy for a long time run under PRP and developing a general and chart-based tool for the problem, making it easier to solve the day-to-day practice and operation of equipment. As a result, a generalized chart was developed to support maintenance decisions through the elaboration of an original isometric table and complemented with a step-by-step methodology to determine the optimum time in which the preventive maintenance activities must be implemented. In most cases, these types of maintenance activities will consider a replacement activity.


---

## Body

## 1. Introduction

Preventive maintenance is primarily intended to avoid high costs incurred due to unplanned maintenance activities, which usually causes higher downtimes, direct costs, and most times catastrophic determining events for business continuity.Because of the need for planning PM activities at the industry level, a direct requirement emerges in consequence: knowing when to perform these activities. Since the PM has considerable importance in safety-critical processes, where the occurrence of a failure has a high impact, most of the developed models are intended to replace the equipment [1–3]. In this sense, most of the policies are based on constant age replacement of equipment, i.e., constant wear time after the equipment entry into operation; this occurs because all devices are under the same condition when they enter as new. Considering that, in the industry, the amount of equipment that only has replacement as a solution for the failure is frequently high, the operators and engineers related to this equipment have to do a lot of exhaustive calculations to find the best policies of replacement. Most of those aforementioned operations, which are also cited in this paper, involve at least calculating nonclosed integrals by Riemann sum.Given the above considerations, this research focuses on analyzing an optimization model with constant age PM on repairable systems and getting a simple way to approach the procedures involved in this optimization process.In the maintenance management context, a repairable system can be defined as “a system that after failing to develop an activity can be restored properly functioning satisfactorily by some method” [4].Depending on the type of maintenance that is given to the equipment, it is possible to find 5 cases [5,6]:(a)
Perfect maintenance: maintenance operation restores the equipment to a status of “as good as new.”(b)
Minimum maintenance: maintenance operation restores the equipment to a status of “as bad as old.”(c)
Imperfect maintenance: maintenance operation restores the equipment to a status of “worse than new but better than old.”(d)
More than perfect maintenance: maintenance operation restores the equipment to a status of “better than new.”(e)
Destructive maintenance: maintenance operation restores the equipment to a status of “worse than old.”Cases (d) and (e) have little presence in reality [6]; almost all correspond to specific events in the equipment life and are highly unlikely to present a series of repairs to such standard.In case (a), the most widely used model that fits precisely the described scenario is the Perfect Renewal Process (PRP) [6,7]; in case (b) a Nonhomogenous Poisson Process (NHPP) [6,7] is used, among others. In the case of (c), various models of failure processes, including variants of generalized geometric processes (for details, see [8]), have been developed; if, besides, the effect of the repair is considered to be limited due to a constant value, then the general repair model is obtained and applied for this particular case [9]; in the same way, the Generalized Renewal Process (GRP) model developed by Kijima and Sumita [10] has been used.The current intense international competition, the need for faster innovation, and the structural organizational changes in response to new information and manufacturing technologies are still surprising [11]. In response to this phenomenon, the increasing importance of effective maintenance planning is fueled by increasing automation of processes, so that a substantial part of the total costs of industrial assets is today directly or indirectly related to maintenance [12].In this sense, maintenance personnel need simple and general methods for the design of maintenance programs and optimization policies [13]. Also, the decision-making tools are needed to give quick but precise answers, simplifying an intricate application of the theory, capable of offering standardized values and ideally easy to read, execute, and interpret.Considering the motivations established above, the main contribution of this research is the elaboration of a practical tool for obtaining fast and reliable information for the decision-making process in maintenance management, by developing a generalized chart-based tool to obtain optimal preventive maintenance times in single-unit systems. Based on the knowledge gathered from the literature review, there are no references related to the development and use of chart-based tools that allow a practical approach to obtain optimal preventive maintenance times.Besides, it is possible to realize that the models developed for equipment replacement can only be used for elements maintained under Perfect Renewal Process modeling, considering that their condition is “as good as new” once they begin its operational functions.The research is structured as follows: Section2 states the problem and its necessity in the industry, Section 3 presents the state of the art related to the PRP modeling and the chosen long-run cost minimization model. Section 4 develops the generalized chart-based proposal solution for the problem, starting with a mathematical analysis of the problem and ends with presenting the results. Finally, Section 5 presents the discussion, potentialities, and conclusions.


## 2. Problem Statement and Motivation

Specialized literature reveals that the industry continues to demand methodologies, approaches, and tools to better understand machine inefficiencies [14]. Specifically, in the maintenance management area, the continuous development of technical systems and the increasing dependence on equipment have generated an increasing interest in effective maintenance activities [15]. In this sense, it becomes necessary not only to determine the maintenance activity to be executed but also to optimize the execution intervals of said activities and correctly model the random behavior of the failures.This approach is targeted at the traditional Perfect Renewal Process (PRP), which is a complement to the existent solution procedures and is also easy to understand when searching and getting the optimal PM time value. Although this stochastic process for the random modeling of failures is elementary, it is the starting point for the elaboration and conceptualization of more sophisticated, but also more complex, models.According to [16], if the purpose of maintenance models is to increase the knowledge of the data, this knowledge should be organized so that the end-user, maintenance personnel, can use it. Current models for getting the optimal time of a PM activity imply exhaustive calculation to get solutions, mainly because of the evaluation of the objective function to be minimized. This means that when it is required to establish the optimization procedure for the equipment, the computation operations of these values can become expensive in terms of time and resources used for this task. Something that also must be considered is the fact that the nature of the problem is completely analytical, requiring an advanced domain of calculus and engineering skills by the user to face and solve the problem.Considering the previous observations, the motivation that guides this research is to find a practical solution to be used quickly, which ideally uses graphical support to find a quick convergence process for an optimal solution.


## 3. State of the Art

### 3.1. Literature Review

Maintenance management has become increasingly relevant within business management, becoming one of its essential functions and also one of the critical aspects within the operational phase of asset management given the current conditions and needs of a highly competitive market and the increasing automation of its processes. In this context, maintenance optimization consists of the development and analysis of mathematical models to optimize maintenance policies [15]; this process is, therefore, essential for the decision-making process in asset management.In this sense, the correct modeling of the failure process is relevant, as well as the different units or components of the system on which the maintenance activities will be carried out. A wide variety of failure process models have been developed and extended over the past three decades addressing different characteristics of failure processes, namely, geometric processes [17], renovation processes [10], virtual age models [18], arithmetic models of intensity reduction (ARI), and arithmetic reduction of age (ARA) [19], among others. To determine the different extensions and generalizations that have been made, Wu et al. [8] elaborate a review of the extensions of the geometric processes, from which it follows that the general repair models are derived from a particular case of the generalized geometric processes.Regarding the relevant investigations carried out in recent years, the research developed by Mullor et al. [20] considers establishing the failure distribution experimentally, disaggregating a specific and independent model for each component of the equipment, and also associating an imperfect maintenance model with each case. In other words, the model developed in said research considers multiple failure modes depending on the different components that integrate the equipment, addressing the implicit assumption that modeling considers a unique degradation process.To address this assumption in turn, [21] develops a multi-dependent system degradation process modeled by a copula function and influenced by the environment, using finite analysis for performance evaluation of life cycle cost. While these studies address single-unit systems, [22] proposes a new failure process model for a series component system. For this problem, a distinction is made between the real system, the mathematical model, and the system made up of virtual components.The previous research works address maintenance activities considering that the failure is self-announcing. The research carried out by Liu et al. [23] considers addressing condition-based maintenance (CBM) policy, and therefore the failure or condition of the equipment can only be detected under the inspection, which is referred to as a soft failure. The problem addressed considers both a system and a sensor in a degraded state, so the traditional maintenance policy without considering this fact may be operating in a suboptimal state.The investigations detailed above are based on the development and application of stochastic processes and often multiple-objective or cost optimization models, with all the complexity that this involves. Based on the classifications used in the literature of the field, the proposed research is framed on the development of a practical tool to obtain optimal maintenance preventive maintenance intervals for single-unit systems and thus a single degradation process, considering perfect maintenance for preventive maintenance activities and minimal maintenance for corrective maintenance actions. Considering the above, said interval is calculated based on the optimization of the expected cost in the long term for said activities.


### 3.2. Maintenance Tools for Decision Making

Related to the problem previously stated, other graphical tools for supporting decision making in the area of operational reliability have been developed.One of the classic management tools for decision making is Pareto diagram. Although this tool is widely used in management, this diagram is not only used to prioritize risk in quality control but also to get Pareto-efficient borders to optimize maintenance costs. In this sense, this tool allows considering the impact of costs from spare parts to downtime and outages [24]. Since the available resources are always scarce and considering the increasing search for integrating different areas within the organization, the search for a Pareto-efficient solution can even be carried out to address a multicriteria objective function [25].The GAMM (Graphical Analysis for Maintenance Management) is a method and maintenance management decision-making tool to visualize and analyze equipment dependability data in a graphical form [26]. One of the advantages of this method is that it supports decision making in the overall maintenance management, through the visualization and graphic analysis of data; for the application of the method, a relatively simple database is enough to generate the graphic support of this tool, consisting in scatter diagrams, bubble charts, and graphs to measure the correlation between each of the equipment reviews and the existing level of equipment reliability before these [27].The Graphical Analysis for Operation Management (GAOM) is a maintenance and operation management tool that is also built from an integrated database but considers three main types of information: maintenance interventions, production, and time. All this information is processed together with values of the key performance indicators for the subsequent development of an integrated and exhaustive graphic analysis of the decision-making process [28].Reviewing the scientific literature, we found that the use of generalized chart-based tools is applied in a wide variety of research fields (see, e.g., [29–32]). This reveals the existing potential and effective utilization of this type of instruments, which can be exploited in the field of maintenance management.


### 3.3. Perfect Renewal Process

The modeling under the Perfect Renewal Process is characterized by considering the perfect operational condition of the equipment (AGAN) after maintenance activity execution. This means that the maintenance activity was fully effective in the rejuvenation process of the equipment such that it returns the equipment to its original condition [6, 7, 10].A complete renovation implies that the time between failures or any unit for which reliability is calculated has no trend over the total time and is independent and identically distributed (i.i.d.) throughout the entire planning horizon.The PRP modeling is the simplest of all, but at the same time it is the one with the greatest limitations in terms of its implementation, not faithfully representing the behavior of most mechanical and electrical equipment. However, PRP modeling is ideal for items that cannot be repaired, considering that after failing there is no way to restore the functionality of these items, and therefore they must be replaced. Thus, whenever a failure occurs, the item is replaced by a new one, describing exactly the operational scenario of PRP. In this model, the failure rate is characterized by the following expression:(1)λtt≥T=ft−T1−Ft−T.Here, ft corresponds to the p.d.f. of the operating times, λT the failure rate at the elapsed time T since the start of the operation, and ti is the i-th operation time.In the case applied to the 2-parameter Weibull distribution, the basic expressions are(2)Ft=1−e−t/αβ,Rt=e−t/αβ,ft=βαtαβ−1e−t/αβ.A graphical representation of how failures occur under this modeling is shown in Figure1, which seeks to represent the independence in the occurrence time of the failures.Figure 1
Failure rate in PRP.Through the use of the inverse function, it is possible to generate a random time value with the following expression:(3)ti=α−lnuβ,where ti represents the time between failures and u a uniformly random variable in the [0, 1] interval.


### 3.4. Optimal Preventive Maintenance in PRP

The Wang and Christer model [33] is the most widely accepted and is based on the minimization of the expected cost in the long run. For this analysis, it is required to consider two possibilities:(a)
The equipment fails before performing PM(b)
The equipment survives until it has performed the maintenance activitySince under modeling PRP, times are identically and independently distributed, and the p.d.f. of the times between failures is the same for all cases, which is also the main difference from NHPP and GRP repairable models. This feature makes the probability of occurrence for both the first and second cases the same over time, so there is a single optimum time in which the PM activity is performed.The occurrence of each of these scenarios implies incurring certain maintenance operations, which have a particular cost and time associated with each case. Namely, in case of no failing before PM (case (b)), the equipment works continuously up to an arbitrary intervaltp, executing only the PM activity immediately after this event. The above can be represented mathematically by establishing that the random variable related to the time of occurrence of the failure tf takes a greater value than the moment of execution of the PM activity tp. Meanwhile, in case (a), the occurrence of the early failure does not allow the preventive maintenance activity to be executed at instant tp, being necessary to carry out repair tasks at the instant of time tf, in which the failure occurs. Mathematically speaking, time becomes a random variable called tf and is represented by the truncated distribution of the times between failures, bounded above by value tp. Figure 2 provides a summary of the previously stated cases. Additionally, the probability associated with each case is clearly shown in Figure 3.Figure 2
Cases representation.Figure 3
p.d.f. of time between failures.Since in case (b) the operating time istp, for case (a) the expected time of operation is given by(4)Ett<tp=∫0tpt⋅ftdtFtp.The model also considers a differentiated cost for each of the two cases, where the cost of case (a) is usually larger than case (b), being an unplanned task. This definesCc presented to the cost of case (a) (corrective cost) and Cp for case (b) (preventive cost).Considering that each of the cases has different lengths, the model suggests minimizing the cost per time rate ([$/t]).For a certaintp, the expected cost of C is given by the following expression:(5)EC=Ftp⋅Cc︸Casea+Rtp⋅Cp︸Caseb.Besides, the expected duration after each maintenance activity is(6)Et=Ftp⋅∫0tpt⋅ftdtftp︸Casea+Rtp⋅tp︸Caseb.(7)Et=∫0tpt⋅ftdt︸Casea+Rt⋅tp︸Caseb.Thus, the expected cost per unit of time can be expressed as(8)EdCttpdt=Ftp⋅Cc+Rtp⋅Cp∫0tpt⋅ftdt+Rtp⋅tp.


## 3.1. Literature Review

Maintenance management has become increasingly relevant within business management, becoming one of its essential functions and also one of the critical aspects within the operational phase of asset management given the current conditions and needs of a highly competitive market and the increasing automation of its processes. In this context, maintenance optimization consists of the development and analysis of mathematical models to optimize maintenance policies [15]; this process is, therefore, essential for the decision-making process in asset management.In this sense, the correct modeling of the failure process is relevant, as well as the different units or components of the system on which the maintenance activities will be carried out. A wide variety of failure process models have been developed and extended over the past three decades addressing different characteristics of failure processes, namely, geometric processes [17], renovation processes [10], virtual age models [18], arithmetic models of intensity reduction (ARI), and arithmetic reduction of age (ARA) [19], among others. To determine the different extensions and generalizations that have been made, Wu et al. [8] elaborate a review of the extensions of the geometric processes, from which it follows that the general repair models are derived from a particular case of the generalized geometric processes.Regarding the relevant investigations carried out in recent years, the research developed by Mullor et al. [20] considers establishing the failure distribution experimentally, disaggregating a specific and independent model for each component of the equipment, and also associating an imperfect maintenance model with each case. In other words, the model developed in said research considers multiple failure modes depending on the different components that integrate the equipment, addressing the implicit assumption that modeling considers a unique degradation process.To address this assumption in turn, [21] develops a multi-dependent system degradation process modeled by a copula function and influenced by the environment, using finite analysis for performance evaluation of life cycle cost. While these studies address single-unit systems, [22] proposes a new failure process model for a series component system. For this problem, a distinction is made between the real system, the mathematical model, and the system made up of virtual components.The previous research works address maintenance activities considering that the failure is self-announcing. The research carried out by Liu et al. [23] considers addressing condition-based maintenance (CBM) policy, and therefore the failure or condition of the equipment can only be detected under the inspection, which is referred to as a soft failure. The problem addressed considers both a system and a sensor in a degraded state, so the traditional maintenance policy without considering this fact may be operating in a suboptimal state.The investigations detailed above are based on the development and application of stochastic processes and often multiple-objective or cost optimization models, with all the complexity that this involves. Based on the classifications used in the literature of the field, the proposed research is framed on the development of a practical tool to obtain optimal maintenance preventive maintenance intervals for single-unit systems and thus a single degradation process, considering perfect maintenance for preventive maintenance activities and minimal maintenance for corrective maintenance actions. Considering the above, said interval is calculated based on the optimization of the expected cost in the long term for said activities.


## 3.2. Maintenance Tools for Decision Making

Related to the problem previously stated, other graphical tools for supporting decision making in the area of operational reliability have been developed.One of the classic management tools for decision making is Pareto diagram. Although this tool is widely used in management, this diagram is not only used to prioritize risk in quality control but also to get Pareto-efficient borders to optimize maintenance costs. In this sense, this tool allows considering the impact of costs from spare parts to downtime and outages [24]. Since the available resources are always scarce and considering the increasing search for integrating different areas within the organization, the search for a Pareto-efficient solution can even be carried out to address a multicriteria objective function [25].The GAMM (Graphical Analysis for Maintenance Management) is a method and maintenance management decision-making tool to visualize and analyze equipment dependability data in a graphical form [26]. One of the advantages of this method is that it supports decision making in the overall maintenance management, through the visualization and graphic analysis of data; for the application of the method, a relatively simple database is enough to generate the graphic support of this tool, consisting in scatter diagrams, bubble charts, and graphs to measure the correlation between each of the equipment reviews and the existing level of equipment reliability before these [27].The Graphical Analysis for Operation Management (GAOM) is a maintenance and operation management tool that is also built from an integrated database but considers three main types of information: maintenance interventions, production, and time. All this information is processed together with values of the key performance indicators for the subsequent development of an integrated and exhaustive graphic analysis of the decision-making process [28].Reviewing the scientific literature, we found that the use of generalized chart-based tools is applied in a wide variety of research fields (see, e.g., [29–32]). This reveals the existing potential and effective utilization of this type of instruments, which can be exploited in the field of maintenance management.


## 3.3. Perfect Renewal Process

The modeling under the Perfect Renewal Process is characterized by considering the perfect operational condition of the equipment (AGAN) after maintenance activity execution. This means that the maintenance activity was fully effective in the rejuvenation process of the equipment such that it returns the equipment to its original condition [6, 7, 10].A complete renovation implies that the time between failures or any unit for which reliability is calculated has no trend over the total time and is independent and identically distributed (i.i.d.) throughout the entire planning horizon.The PRP modeling is the simplest of all, but at the same time it is the one with the greatest limitations in terms of its implementation, not faithfully representing the behavior of most mechanical and electrical equipment. However, PRP modeling is ideal for items that cannot be repaired, considering that after failing there is no way to restore the functionality of these items, and therefore they must be replaced. Thus, whenever a failure occurs, the item is replaced by a new one, describing exactly the operational scenario of PRP. In this model, the failure rate is characterized by the following expression:(1)λtt≥T=ft−T1−Ft−T.Here, ft corresponds to the p.d.f. of the operating times, λT the failure rate at the elapsed time T since the start of the operation, and ti is the i-th operation time.In the case applied to the 2-parameter Weibull distribution, the basic expressions are(2)Ft=1−e−t/αβ,Rt=e−t/αβ,ft=βαtαβ−1e−t/αβ.A graphical representation of how failures occur under this modeling is shown in Figure1, which seeks to represent the independence in the occurrence time of the failures.Figure 1
Failure rate in PRP.Through the use of the inverse function, it is possible to generate a random time value with the following expression:(3)ti=α−lnuβ,where ti represents the time between failures and u a uniformly random variable in the [0, 1] interval.


## 3.4. Optimal Preventive Maintenance in PRP

The Wang and Christer model [33] is the most widely accepted and is based on the minimization of the expected cost in the long run. For this analysis, it is required to consider two possibilities:(a)
The equipment fails before performing PM(b)
The equipment survives until it has performed the maintenance activitySince under modeling PRP, times are identically and independently distributed, and the p.d.f. of the times between failures is the same for all cases, which is also the main difference from NHPP and GRP repairable models. This feature makes the probability of occurrence for both the first and second cases the same over time, so there is a single optimum time in which the PM activity is performed.The occurrence of each of these scenarios implies incurring certain maintenance operations, which have a particular cost and time associated with each case. Namely, in case of no failing before PM (case (b)), the equipment works continuously up to an arbitrary intervaltp, executing only the PM activity immediately after this event. The above can be represented mathematically by establishing that the random variable related to the time of occurrence of the failure tf takes a greater value than the moment of execution of the PM activity tp. Meanwhile, in case (a), the occurrence of the early failure does not allow the preventive maintenance activity to be executed at instant tp, being necessary to carry out repair tasks at the instant of time tf, in which the failure occurs. Mathematically speaking, time becomes a random variable called tf and is represented by the truncated distribution of the times between failures, bounded above by value tp. Figure 2 provides a summary of the previously stated cases. Additionally, the probability associated with each case is clearly shown in Figure 3.Figure 2
Cases representation.Figure 3
p.d.f. of time between failures.Since in case (b) the operating time istp, for case (a) the expected time of operation is given by(4)Ett<tp=∫0tpt⋅ftdtFtp.The model also considers a differentiated cost for each of the two cases, where the cost of case (a) is usually larger than case (b), being an unplanned task. This definesCc presented to the cost of case (a) (corrective cost) and Cp for case (b) (preventive cost).Considering that each of the cases has different lengths, the model suggests minimizing the cost per time rate ([$/t]).For a certaintp, the expected cost of C is given by the following expression:(5)EC=Ftp⋅Cc︸Casea+Rtp⋅Cp︸Caseb.Besides, the expected duration after each maintenance activity is(6)Et=Ftp⋅∫0tpt⋅ftdtftp︸Casea+Rtp⋅tp︸Caseb.(7)Et=∫0tpt⋅ftdt︸Casea+Rt⋅tp︸Caseb.Thus, the expected cost per unit of time can be expressed as(8)EdCttpdt=Ftp⋅Cc+Rtp⋅Cp∫0tpt⋅ftdt+Rtp⋅tp.


## 4. Research Proposal

To develop the methodology, firstly it is required to take the model from the previous point to tangible terms and simplify the number of parameters included. In the case described by two-parameter Weibull distribution, the expected cost per unit of time corresponds to(9)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cp∫0tpt⋅β/αt/αβ−1e−t/αβdt︸k+e−tp/αβ⋅tp.Here, the k variable corresponds to(10)k=∫0tpt⋅βαtαβ−1e−t/αβdt.By changing the variables for the definedk part,(11)tβαβ=u,αβu1/β=t,αβu1/ββudu=dt.Replacing the integration limits,(12)t=tp⟶u=tpαβ,t=0⟶u=0,and considering that t/α=u1/β,(13)k=∫0tp/αβαβu1/β⋅βαu1−1/βe−u⋅αβu1/ββudu,k=∫0tp/αβαβu1/β⋅βαu1−1/βe−u⋅αβu1/ββudu=∫0tp/αβαu1/β⋅βαuu−1/βe−u⋅αu1/ββudu.Thus,k satisfies that(14)k=α∫0tp/αβu1/βe−udu=α∫0tp/αβu1/β+1−1e−udu︸γ1/β+1,tp/αβ,k=α⋅γ1/β+1,tpαβ,where the γa,x function represents the lower incomplete gamma function:(15)γa,x=∫0xua−1du.This leads to formulating expression (8) as follows:(16)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.Then, the objective is to minimize the long-run expected cost; i.e.,(17)Mintp∈ℝ+:1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.Before continuing, it is important to consider that the value oftp is determined by Cc, Cp, β, and α, and it is possible to simplify the expression stated in (17) by using the rate between Cc and Cp and using a single value for α that only impacts the size of tp. It is clear that α and tp have a direct proportionality. In other words, for the ratio between costs for each different β, there is a tp which is directly proportional to α. Considering that it is possible to formulate a change of variable(18)tp′=tpα,this further simplifies the denominator, since it represents the expected duration.At this point, the rate(19)p=CpCc,is defined, with the purpose of diminishing the variables involved in the objective function stated in (17). Finally, the function depends only ontp′ and p, stated as follows:(20)Mintp∈ℝ+:1−e−tp′β⋅Cc+e−tp′β⋅pCpα⋅γ1+1/β,tp′β+e−tp′β⋅tp′.It is also important to note that, using the lower incomplete gamma function, it is possible to rewrite it as(21)γa,x=Pa,x⋅Γa.Here, Pa,x is the cumulative distribution function of a gamma distribution probability with a scale parameter equal to 1 and a shape parameter equal to a. This change greatly improves the possibility of reaching a better solution, given that most calculation software has already defined the function approximations for the probability distribution function of gamma and the gamma function itself. Then, the equivalence between the variables is(22)a=1+1β,x=tp′β.Then, by changing the corresponding variables(23)γ1+1β,tp′β=P1+1β,tp′β⋅Γ1+1β,and applying all the mentioned changes, expression (17) is equal to(24)Mintp′∈ℝ+:Ccα⋅1−e−tp′β⋅1−pP1+1/β,tp′β⋅Γ1+1/β+e−tp′β⋅tp′,which is equivalent to solving the following expression:(25)Mintp′∈ℝ+:1−e−tp′β⋅1−pP1+1/β,tp′β⋅Γ1+1/β+e−tp′β⋅tp′.


### 4.1. Constant Failure Rate Analysis

First, an analysis is made for the value ofβ=1, with a value p=0, so the preventive cost is zero, an extreme case value applied to expression (8):(26)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.By simplifying(27)α=1,p=CpCc=0,β=1,and replacing the values(28)EdCttpdt=Cc1−e−tpγ2,tp+e−tp⋅tp,evidently the expected cost minimization for a given time value does not depend on Cc; this is a scale parameter for the problem. Then, the problem becomes(29)Mintp∈ℝ+:1−e−tpγ2,tp+e−tp⋅tp.However,(30)γ2,t+e−t⋅t=∫0ts2−1e−sds︸w+e−t⋅t.Taking the indefinite integral of thew part(31)w=∫s2−1e−sds,it is possible to develop the integral by parts∫udv=uv−∫vdu, making the changes of variables(32)u=s,dv=e−sds,du=ds,v=−e−s.This leads to(33)w=−e−ss+∫e−sds=−e−ss−e−s+C,in the definite form(34)∫0ts2−1e−sds=−e−ss−et0=−e−tt+1+1.So(35)γ2,t+e−t⋅t=∫0ts2−1e−sds+e−1⋅t=−e−tt+1+1+e−t⋅t=1−e−t.Therefore, expression (29) is equal to(36)1−e−tpγ2,tp+e−tp×tp=1−e−tp1−e−tp=1∀tp>0.The latter demonstrates and verifies that it is indifferent to the time when PM is performed: the cost will be the same. With this result, it is possible to establish that if the PM cost is not null but less than the cost of corrective maintenance, it will never be desirable to carry out a PM activity, so it is better to wait for the occurrence of the failure since it will decrease the cost per operating time. This because it has just been proven that a PM cost of zero is indifferent to when PM is executed.


### 4.2. Nonconstant Failure Rate Analysis

Again, we start from expression (8) as:(37)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.We consider a special case, where(38)α=1,p=CpCc=0.Then, we replace the values in expression (8):(39)EdCttpdt=Cc1−e−tpβγ1+1/β,tpβ+e−tβ⋅tp.Based on the above, the value oftp that satisfies the slightest expected cost per unit of time is the same as the solution of the following expression:(40)Mintp∈ℝ+:1−e−tpβγ1+1/β,tpβ+e−tpβ⋅tp.To mathematically analyze this problem is quite complex, mainly due to an iterative integration problem. In this sense, it is interesting to note that if the shape parameter of the Weibull distribution tends to zero, then(41)limβ⟶0γ1+1β,tpβ=0.Therefore, the minimization of the objective function expression becomes trivial, considering an expected cost that tends to a null value. This means that maintenance activities are avoided since, for β∈0,1, they are unproductive. This situation could be resolved by a search algorithm in space β>1, which estimates that the best way to solve the situation is to develop an algorithm to find the optimal time directly for the original function and obtain a graphical tool with the results of this algorithm. The development of this tool is presented in the next section.


### 4.3. Development of the Generalized Charts

Considering the defined requirements of the problem in the above sections, it has been justified that a graphical tool is really necessary. The problem has been simplified, and it is now possible to generalize the problem consideringβ, p, and α. This amount of variables tells us that it could be possible to plot all the solutions of the problems on a chart, which is the challenge of this section.


#### 4.3.1. Search Algorithm

With regard to the function already defined, it must be remembered that the parameters that alter the problem areβ and p, given that α is a scale parameter that is extracted and generalized with value 1. In other words, without loss of generality, it is possible to use for this purpose the expression in (25) considering that tp=tp′. Once these optimal values for this standard tp value are obtained, then it is enough to only scale the value from the estimation of the shape parameter α to obtain the practical optimal PM time value. Establishing the previous observations, we have developed an algorithm for finding optimal PM time for different values of β and p.The time exploring ranges fromt=0 to t=10−5α. With α=1, the initial value is t=10−5 time units. The maximum scan time is given by the value of t that contains the 1−δ cases for the probability distribution function, t=F−11−δ. Then, δ would be the percentage of the highest values of t that are not being considered for the analysis. Thus, temax value as the maximum exploring time remains as(42)temax=F−11−δ=α−lnδβ.This expression is represented in Figure4.Figure 4
Explored space in the c.d.f. by the algorithm.As expected,δ should have a small value for the sample to be valid. In this research case, 99.9% of the cases have been explored, which implies a value of δ=0.001. The steps in the value of t are also defined, which defines the resolution of the search. Also, certain steps numbers, n1 and n2, are defined, each one with a different space search as shown in Figure 5.Figure 5
Search algorithm.As shown in the previous illustration, the first search space is defined at interval0,F−11−δ. This space is explored by a number of discrete steps n1, which is given by(43)n1=F−11−δr1,where r1 corresponds to the range for each consultation of a new value for the objective function. The second search space has new n2 and r2, corresponding to the number of values to be explored and the of value the range between each observation. A value of n1=n2=103 was used for this study, with a final resolution of a rate close to 10-6.It is important to note that this algorithm is tunable in terms of the size of sample spaces and the number of assessments in space in order to reduce the total number of assessments and for overall better resolution. For this publication, this is not very significant and could be studied in another case where the objective is to optimize the optimization algorithm.Once an optimal value is found in the first search space andtopi value is obtained, we proceed to define the second search space corresponding to segment topi−r1,topi+r1. This range is again explored in n2 equidistant points. Within this range, it is shown that the lowest value of the objective function value is defined as the smallest topt found in the second search space as the optimal value for preventive intervention to the original α value, which for simplification has been α = 1. This topt has been denominated standard optimal time to PM.


#### 4.3.2. Results of Algorithm for PRP

The previous algorithm has been run for different values ofβ and p as indicated in 4.3.1. A 3D-chart representation of several solutions has been obtained and presented in Figure 6.Figure 6
Optimal preventive time for PRP modeling.The chart in Figure6 has been modified, so that areas where the optimal value of PM time gives inconsistent results, such as an undefined value or a value greater than the MTBF, are not considered. Besides, the grid used for the graph is reduced with respect to the resolution of the generalized chart to improve the visualization. Although this figure represents a first approach prior to obtaining generalized charts, important conclusions were obtained from the result:1
There is a clearly defined area where it is convenient to perform PM on modeling PRP.2
There is no time for PM for values ofβ < 1. It is not convenient to perform PM, independently of the ratio between the cost and the cost preventive correction, p.Under the criterion, preventive activity is not performed within 0.1% of the highest times; the results have shown that, in in 99.91% of cases, the optimal values for a givenβ and p are within the search range, previously established by the δ variable. The results are then consistent with the assumptions used.In Figure7, the boundary zone to run optimal PM is also presented. The entire area under the dashed curve corresponds to the pair of values for p and β, in which it is convenient to implement the PM policy in a fixed time, which is considered throughout this work. The area above the curve corresponds to the possibilities where it is not convenient to have a policy of fixed-time PM.Figure 7
Limit possibilities for preventive maintenance under PRP modeling.Finally, a generalized chart has been developed as a result of this research, which is shown schematically in Figure8 and presented in detail in the appendix. Since the most common values for the beta parameter are β∈1,2, the generalized charts for this range are incorporated in detail in Appendix A, using a 0.01 step for each adjustment. For general purposes, generalized charts for β∈1,5 are also incorporated in Appendix B considering a step of 0.05 for each β value.Figure 8
Chart representation for optimal values searching.This chart contains a generalization of the optimization for different parameter values. In other words, to obtain the optimal time of replacement for a 2-parameter Weibull distributed lifetime of an item, first, thep parameter is needed, i.e., the relation between PM cost and corrective maintenance cost. Then, it is necessary to find the corresponding curve for the β value of the distribution. Once the intersection is located, it is necessary to scroll horizontally to the left to find the corresponding topt value. This t value is generalized considering α=1, and then the tp optimal value will be tp=α⋅topt. If there is no intersection between the p parameter and theβ corresponding curve, it is not convenient to perform a PM activity.


## 4.1. Constant Failure Rate Analysis

First, an analysis is made for the value ofβ=1, with a value p=0, so the preventive cost is zero, an extreme case value applied to expression (8):(26)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.By simplifying(27)α=1,p=CpCc=0,β=1,and replacing the values(28)EdCttpdt=Cc1−e−tpγ2,tp+e−tp⋅tp,evidently the expected cost minimization for a given time value does not depend on Cc; this is a scale parameter for the problem. Then, the problem becomes(29)Mintp∈ℝ+:1−e−tpγ2,tp+e−tp⋅tp.However,(30)γ2,t+e−t⋅t=∫0ts2−1e−sds︸w+e−t⋅t.Taking the indefinite integral of thew part(31)w=∫s2−1e−sds,it is possible to develop the integral by parts∫udv=uv−∫vdu, making the changes of variables(32)u=s,dv=e−sds,du=ds,v=−e−s.This leads to(33)w=−e−ss+∫e−sds=−e−ss−e−s+C,in the definite form(34)∫0ts2−1e−sds=−e−ss−et0=−e−tt+1+1.So(35)γ2,t+e−t⋅t=∫0ts2−1e−sds+e−1⋅t=−e−tt+1+1+e−t⋅t=1−e−t.Therefore, expression (29) is equal to(36)1−e−tpγ2,tp+e−tp×tp=1−e−tp1−e−tp=1∀tp>0.The latter demonstrates and verifies that it is indifferent to the time when PM is performed: the cost will be the same. With this result, it is possible to establish that if the PM cost is not null but less than the cost of corrective maintenance, it will never be desirable to carry out a PM activity, so it is better to wait for the occurrence of the failure since it will decrease the cost per operating time. This because it has just been proven that a PM cost of zero is indifferent to when PM is executed.


## 4.2. Nonconstant Failure Rate Analysis

Again, we start from expression (8) as:(37)EdCttpdt=1−e−tp/αβ⋅Cc+e−tp/αβ⋅Cpα⋅γ1+1/β,tp/αβ+e−tp/αβ⋅tp.We consider a special case, where(38)α=1,p=CpCc=0.Then, we replace the values in expression (8):(39)EdCttpdt=Cc1−e−tpβγ1+1/β,tpβ+e−tβ⋅tp.Based on the above, the value oftp that satisfies the slightest expected cost per unit of time is the same as the solution of the following expression:(40)Mintp∈ℝ+:1−e−tpβγ1+1/β,tpβ+e−tpβ⋅tp.To mathematically analyze this problem is quite complex, mainly due to an iterative integration problem. In this sense, it is interesting to note that if the shape parameter of the Weibull distribution tends to zero, then(41)limβ⟶0γ1+1β,tpβ=0.Therefore, the minimization of the objective function expression becomes trivial, considering an expected cost that tends to a null value. This means that maintenance activities are avoided since, for β∈0,1, they are unproductive. This situation could be resolved by a search algorithm in space β>1, which estimates that the best way to solve the situation is to develop an algorithm to find the optimal time directly for the original function and obtain a graphical tool with the results of this algorithm. The development of this tool is presented in the next section.


## 4.3. Development of the Generalized Charts

Considering the defined requirements of the problem in the above sections, it has been justified that a graphical tool is really necessary. The problem has been simplified, and it is now possible to generalize the problem consideringβ, p, and α. This amount of variables tells us that it could be possible to plot all the solutions of the problems on a chart, which is the challenge of this section.


### 4.3.1. Search Algorithm

With regard to the function already defined, it must be remembered that the parameters that alter the problem areβ and p, given that α is a scale parameter that is extracted and generalized with value 1. In other words, without loss of generality, it is possible to use for this purpose the expression in (25) considering that tp=tp′. Once these optimal values for this standard tp value are obtained, then it is enough to only scale the value from the estimation of the shape parameter α to obtain the practical optimal PM time value. Establishing the previous observations, we have developed an algorithm for finding optimal PM time for different values of β and p.The time exploring ranges fromt=0 to t=10−5α. With α=1, the initial value is t=10−5 time units. The maximum scan time is given by the value of t that contains the 1−δ cases for the probability distribution function, t=F−11−δ. Then, δ would be the percentage of the highest values of t that are not being considered for the analysis. Thus, temax value as the maximum exploring time remains as(42)temax=F−11−δ=α−lnδβ.This expression is represented in Figure4.Figure 4
Explored space in the c.d.f. by the algorithm.As expected,δ should have a small value for the sample to be valid. In this research case, 99.9% of the cases have been explored, which implies a value of δ=0.001. The steps in the value of t are also defined, which defines the resolution of the search. Also, certain steps numbers, n1 and n2, are defined, each one with a different space search as shown in Figure 5.Figure 5
Search algorithm.As shown in the previous illustration, the first search space is defined at interval0,F−11−δ. This space is explored by a number of discrete steps n1, which is given by(43)n1=F−11−δr1,where r1 corresponds to the range for each consultation of a new value for the objective function. The second search space has new n2 and r2, corresponding to the number of values to be explored and the of value the range between each observation. A value of n1=n2=103 was used for this study, with a final resolution of a rate close to 10-6.It is important to note that this algorithm is tunable in terms of the size of sample spaces and the number of assessments in space in order to reduce the total number of assessments and for overall better resolution. For this publication, this is not very significant and could be studied in another case where the objective is to optimize the optimization algorithm.Once an optimal value is found in the first search space andtopi value is obtained, we proceed to define the second search space corresponding to segment topi−r1,topi+r1. This range is again explored in n2 equidistant points. Within this range, it is shown that the lowest value of the objective function value is defined as the smallest topt found in the second search space as the optimal value for preventive intervention to the original α value, which for simplification has been α = 1. This topt has been denominated standard optimal time to PM.


### 4.3.2. Results of Algorithm for PRP

The previous algorithm has been run for different values ofβ and p as indicated in 4.3.1. A 3D-chart representation of several solutions has been obtained and presented in Figure 6.Figure 6
Optimal preventive time for PRP modeling.The chart in Figure6 has been modified, so that areas where the optimal value of PM time gives inconsistent results, such as an undefined value or a value greater than the MTBF, are not considered. Besides, the grid used for the graph is reduced with respect to the resolution of the generalized chart to improve the visualization. Although this figure represents a first approach prior to obtaining generalized charts, important conclusions were obtained from the result:1
There is a clearly defined area where it is convenient to perform PM on modeling PRP.2
There is no time for PM for values ofβ < 1. It is not convenient to perform PM, independently of the ratio between the cost and the cost preventive correction, p.Under the criterion, preventive activity is not performed within 0.1% of the highest times; the results have shown that, in in 99.91% of cases, the optimal values for a givenβ and p are within the search range, previously established by the δ variable. The results are then consistent with the assumptions used.In Figure7, the boundary zone to run optimal PM is also presented. The entire area under the dashed curve corresponds to the pair of values for p and β, in which it is convenient to implement the PM policy in a fixed time, which is considered throughout this work. The area above the curve corresponds to the possibilities where it is not convenient to have a policy of fixed-time PM.Figure 7
Limit possibilities for preventive maintenance under PRP modeling.Finally, a generalized chart has been developed as a result of this research, which is shown schematically in Figure8 and presented in detail in the appendix. Since the most common values for the beta parameter are β∈1,2, the generalized charts for this range are incorporated in detail in Appendix A, using a 0.01 step for each adjustment. For general purposes, generalized charts for β∈1,5 are also incorporated in Appendix B considering a step of 0.05 for each β value.Figure 8
Chart representation for optimal values searching.This chart contains a generalization of the optimization for different parameter values. In other words, to obtain the optimal time of replacement for a 2-parameter Weibull distributed lifetime of an item, first, thep parameter is needed, i.e., the relation between PM cost and corrective maintenance cost. Then, it is necessary to find the corresponding curve for the β value of the distribution. Once the intersection is located, it is necessary to scroll horizontally to the left to find the corresponding topt value. This t value is generalized considering α=1, and then the tp optimal value will be tp=α⋅topt. If there is no intersection between the p parameter and theβ corresponding curve, it is not convenient to perform a PM activity.


## 4.3.1. Search Algorithm

With regard to the function already defined, it must be remembered that the parameters that alter the problem areβ and p, given that α is a scale parameter that is extracted and generalized with value 1. In other words, without loss of generality, it is possible to use for this purpose the expression in (25) considering that tp=tp′. Once these optimal values for this standard tp value are obtained, then it is enough to only scale the value from the estimation of the shape parameter α to obtain the practical optimal PM time value. Establishing the previous observations, we have developed an algorithm for finding optimal PM time for different values of β and p.The time exploring ranges fromt=0 to t=10−5α. With α=1, the initial value is t=10−5 time units. The maximum scan time is given by the value of t that contains the 1−δ cases for the probability distribution function, t=F−11−δ. Then, δ would be the percentage of the highest values of t that are not being considered for the analysis. Thus, temax value as the maximum exploring time remains as(42)temax=F−11−δ=α−lnδβ.This expression is represented in Figure4.Figure 4
Explored space in the c.d.f. by the algorithm.As expected,δ should have a small value for the sample to be valid. In this research case, 99.9% of the cases have been explored, which implies a value of δ=0.001. The steps in the value of t are also defined, which defines the resolution of the search. Also, certain steps numbers, n1 and n2, are defined, each one with a different space search as shown in Figure 5.Figure 5
Search algorithm.As shown in the previous illustration, the first search space is defined at interval0,F−11−δ. This space is explored by a number of discrete steps n1, which is given by(43)n1=F−11−δr1,where r1 corresponds to the range for each consultation of a new value for the objective function. The second search space has new n2 and r2, corresponding to the number of values to be explored and the of value the range between each observation. A value of n1=n2=103 was used for this study, with a final resolution of a rate close to 10-6.It is important to note that this algorithm is tunable in terms of the size of sample spaces and the number of assessments in space in order to reduce the total number of assessments and for overall better resolution. For this publication, this is not very significant and could be studied in another case where the objective is to optimize the optimization algorithm.Once an optimal value is found in the first search space andtopi value is obtained, we proceed to define the second search space corresponding to segment topi−r1,topi+r1. This range is again explored in n2 equidistant points. Within this range, it is shown that the lowest value of the objective function value is defined as the smallest topt found in the second search space as the optimal value for preventive intervention to the original α value, which for simplification has been α = 1. This topt has been denominated standard optimal time to PM.


## 4.3.2. Results of Algorithm for PRP

The previous algorithm has been run for different values ofβ and p as indicated in 4.3.1. A 3D-chart representation of several solutions has been obtained and presented in Figure 6.Figure 6
Optimal preventive time for PRP modeling.The chart in Figure6 has been modified, so that areas where the optimal value of PM time gives inconsistent results, such as an undefined value or a value greater than the MTBF, are not considered. Besides, the grid used for the graph is reduced with respect to the resolution of the generalized chart to improve the visualization. Although this figure represents a first approach prior to obtaining generalized charts, important conclusions were obtained from the result:1
There is a clearly defined area where it is convenient to perform PM on modeling PRP.2
There is no time for PM for values ofβ < 1. It is not convenient to perform PM, independently of the ratio between the cost and the cost preventive correction, p.Under the criterion, preventive activity is not performed within 0.1% of the highest times; the results have shown that, in in 99.91% of cases, the optimal values for a givenβ and p are within the search range, previously established by the δ variable. The results are then consistent with the assumptions used.In Figure7, the boundary zone to run optimal PM is also presented. The entire area under the dashed curve corresponds to the pair of values for p and β, in which it is convenient to implement the PM policy in a fixed time, which is considered throughout this work. The area above the curve corresponds to the possibilities where it is not convenient to have a policy of fixed-time PM.Figure 7
Limit possibilities for preventive maintenance under PRP modeling.Finally, a generalized chart has been developed as a result of this research, which is shown schematically in Figure8 and presented in detail in the appendix. Since the most common values for the beta parameter are β∈1,2, the generalized charts for this range are incorporated in detail in Appendix A, using a 0.01 step for each adjustment. For general purposes, generalized charts for β∈1,5 are also incorporated in Appendix B considering a step of 0.05 for each β value.Figure 8
Chart representation for optimal values searching.This chart contains a generalization of the optimization for different parameter values. In other words, to obtain the optimal time of replacement for a 2-parameter Weibull distributed lifetime of an item, first, thep parameter is needed, i.e., the relation between PM cost and corrective maintenance cost. Then, it is necessary to find the corresponding curve for the β value of the distribution. Once the intersection is located, it is necessary to scroll horizontally to the left to find the corresponding topt value. This t value is generalized considering α=1, and then the tp optimal value will be tp=α⋅topt. If there is no intersection between the p parameter and theβ corresponding curve, it is not convenient to perform a PM activity.


## 5. Case Study

### 5.1. Description of the Process

The case study and application of the purposed methodology are focused on one of the slurry pumps (P1) of a mining process, which is essential for the continuity of the process (critical) and is exposed to a high-stress job, especially by the transfer characteristics of highly abrasive and corrosive flow, high solids content, and strong acidity.A slurry pump is a rugged heavy-duty pump intended for aggressive or abrasive slurry solutions typically found in the mining industry with particles of various sizes. These pumps increase the pressure of a liquid and solid particle mixture (commonly called slurry) through centrifugal force (a rotating impeller) and convert electrical energy into slurry potential and kinetic energy. In Figure9, a general scheme of the pump is presented.Figure 9
General scheme of slurry pumps.


### 5.2. The Proposed Methodology

In the following are presented the steps of the proposed methodology, as well as the historical collection of failure data to obtain the required parameters, applied in the case study mentioned. (Table1)Step 1
Trend testSlurry pump P1 has registered the following 30 times between failures, measured in hours and presented in order of occurrence.Table 1
Times between failures for slurry pump P1 in hours.
No. of failuresTime between failures (h)No. of failuresTime between failures (h)1434.0916375.052226.0417496.483266.6718649.414681.3319540.6851127.8520491.596634.6921947.317474.8422718.78838.7223953.63931.8924182.210711.5225391.5511726.926327.7112574.1527986.04131043.5428352.5314336.5429631.8615771.2330680.52First of all, as the literature [34] proposes, it is required to evaluate a trend test; the most commonly used is the Laplace Trend Test [34, 35], which gives a test statistic z that has to be compared with the standard normal distribution because it fits approximately to n standard normal random variable. The z value is calculated by(44)z=∑i=1nti/n−to/2to1/12n.Here, to is the total elapsed time, n is the number of times, and ti is the total elapsed time to element i.And the nullH0 and alternative Ha hypotheses consider the following:H0: no trendz=0.Ha:  deteriorating trendz>0.Ha: improving trendz<0.Evaluating the expression gives a value ofz=−0.0953. From the standard normal tables with a significance level of 0.10, the critical value is equal to 1.645. If the z value satisfies −1.645<z<1.645, then we would fail to reject the hypothesis of no trend. Since z=−0.0953, then is not possible to say that a trend exists, so PRP modeling is applicable.Step 2
Distribution fittingOnce PRP modeling is properly justified it is necessary to adjust the parameters for the 2-Weibull distributions for this modeling. The most traditional method is to estimate the parameters maximizing the likelihood. The likelihood function is expressed through the joint probability as(45)Pxiinxi,xi+dx,∀i∈1,…n=∏i=1nfxi;θ,Lθ=∏i=1nfxi;θ,where θ corresponds to the vector of parameters of the distribution related to ft and xi corresponds to the i-th realization of the sample. As the desired maximum is the likelihood between data and a p.d.f.: ftθ, the values of the vector θ have to be adjusted with the objective function to meet that maximum. Conceptually, it is to find the parameters that adjust better to a sample xi,…,xn so that the likelihood of the values provided in a random sampling is maximum: Pxiinxi,xi+dx,∀i∈1,…,n.Then, the parametersθ for this case α,β that provide maximum likelihood are those that meet [36, 37](46)α^=∑i=1ntiβ^n1/β^.(47)∑i=1ntiβ^lnti∑i=1ntiβ^−1β^=1n∑i=1nlnti.Here, ti is the observed time between successive failures and n is the total number of failures.Solving the above equations, the estimators areα^=632.04 and β^=1.99.Step 3
Calculating the cost related parameterIn relation to the cost for the historical failure data, the costs of the corrective interventions have been collected in Table2, expressed in (USD).Table 2
Corrective maintenance cost for each failure in USD.
No. of failuresCost of the failure (USD)No. of failuresCost of the failure (USD)11132.7116847.82972.75171049.033853.5818474.864789.41191740.395836.97201041.646994.3921898.171274.3322594.438528231539.229847.06241496.3101140.66251375.55111010.9526748.0412468.0127877.48131045.2928793.51141493.8529430.25151320.56301375.03The average corrective cost is 999.67 (USD), and it is known that the preventive activity has a cost of 160 (USD), so the value ofp is given by(48)p=c0USD999.67USD≈0.16.Step 4
Chart useThis step is composed of substeps. In this stage, the proposed chart in the appendix is needed. The substeps are the following:(1)
With the value ofp, go to the abscissae axis, in this case 0.16.(2)
Find the correspondingβ value: In this case, we consider 1.99 ≈ 2. Follow the corresponding line for the value of p in vertical sense up to the intersection with the corresponding β curve.(3)
Once the intersection is found, move horizontally to the left and find the standardtopt value.(4)
Calculate the optimaltp time with the product between the standard topt value and α.The steps are schematically represented in Figures10–12.Figure 10
Step 1 for use of the proposed chart.Figure 11
Step 2 for use of the proposed chart.Figure 12
Step 3 for use of the proposed chart.Following Step 4, we find that the optimaltp value when performing the preventive activity is at tp=0.443⋅α=279.99h. Then, the PM policy is ready to be implemented. Notice that, for other values of β, an interpolation could be useful between the curves.


## 5.1. Description of the Process

The case study and application of the purposed methodology are focused on one of the slurry pumps (P1) of a mining process, which is essential for the continuity of the process (critical) and is exposed to a high-stress job, especially by the transfer characteristics of highly abrasive and corrosive flow, high solids content, and strong acidity.A slurry pump is a rugged heavy-duty pump intended for aggressive or abrasive slurry solutions typically found in the mining industry with particles of various sizes. These pumps increase the pressure of a liquid and solid particle mixture (commonly called slurry) through centrifugal force (a rotating impeller) and convert electrical energy into slurry potential and kinetic energy. In Figure9, a general scheme of the pump is presented.Figure 9
General scheme of slurry pumps.


## 5.2. The Proposed Methodology

In the following are presented the steps of the proposed methodology, as well as the historical collection of failure data to obtain the required parameters, applied in the case study mentioned. (Table1)Step 1
Trend testSlurry pump P1 has registered the following 30 times between failures, measured in hours and presented in order of occurrence.Table 1
Times between failures for slurry pump P1 in hours.
No. of failuresTime between failures (h)No. of failuresTime between failures (h)1434.0916375.052226.0417496.483266.6718649.414681.3319540.6851127.8520491.596634.6921947.317474.8422718.78838.7223953.63931.8924182.210711.5225391.5511726.926327.7112574.1527986.04131043.5428352.5314336.5429631.8615771.2330680.52First of all, as the literature [34] proposes, it is required to evaluate a trend test; the most commonly used is the Laplace Trend Test [34, 35], which gives a test statistic z that has to be compared with the standard normal distribution because it fits approximately to n standard normal random variable. The z value is calculated by(44)z=∑i=1nti/n−to/2to1/12n.Here, to is the total elapsed time, n is the number of times, and ti is the total elapsed time to element i.And the nullH0 and alternative Ha hypotheses consider the following:H0: no trendz=0.Ha:  deteriorating trendz>0.Ha: improving trendz<0.Evaluating the expression gives a value ofz=−0.0953. From the standard normal tables with a significance level of 0.10, the critical value is equal to 1.645. If the z value satisfies −1.645<z<1.645, then we would fail to reject the hypothesis of no trend. Since z=−0.0953, then is not possible to say that a trend exists, so PRP modeling is applicable.Step 2
Distribution fittingOnce PRP modeling is properly justified it is necessary to adjust the parameters for the 2-Weibull distributions for this modeling. The most traditional method is to estimate the parameters maximizing the likelihood. The likelihood function is expressed through the joint probability as(45)Pxiinxi,xi+dx,∀i∈1,…n=∏i=1nfxi;θ,Lθ=∏i=1nfxi;θ,where θ corresponds to the vector of parameters of the distribution related to ft and xi corresponds to the i-th realization of the sample. As the desired maximum is the likelihood between data and a p.d.f.: ftθ, the values of the vector θ have to be adjusted with the objective function to meet that maximum. Conceptually, it is to find the parameters that adjust better to a sample xi,…,xn so that the likelihood of the values provided in a random sampling is maximum: Pxiinxi,xi+dx,∀i∈1,…,n.Then, the parametersθ for this case α,β that provide maximum likelihood are those that meet [36, 37](46)α^=∑i=1ntiβ^n1/β^.(47)∑i=1ntiβ^lnti∑i=1ntiβ^−1β^=1n∑i=1nlnti.Here, ti is the observed time between successive failures and n is the total number of failures.Solving the above equations, the estimators areα^=632.04 and β^=1.99.Step 3
Calculating the cost related parameterIn relation to the cost for the historical failure data, the costs of the corrective interventions have been collected in Table2, expressed in (USD).Table 2
Corrective maintenance cost for each failure in USD.
No. of failuresCost of the failure (USD)No. of failuresCost of the failure (USD)11132.7116847.82972.75171049.033853.5818474.864789.41191740.395836.97201041.646994.3921898.171274.3322594.438528231539.229847.06241496.3101140.66251375.55111010.9526748.0412468.0127877.48131045.2928793.51141493.8529430.25151320.56301375.03The average corrective cost is 999.67 (USD), and it is known that the preventive activity has a cost of 160 (USD), so the value ofp is given by(48)p=c0USD999.67USD≈0.16.Step 4
Chart useThis step is composed of substeps. In this stage, the proposed chart in the appendix is needed. The substeps are the following:(1)
With the value ofp, go to the abscissae axis, in this case 0.16.(2)
Find the correspondingβ value: In this case, we consider 1.99 ≈ 2. Follow the corresponding line for the value of p in vertical sense up to the intersection with the corresponding β curve.(3)
Once the intersection is found, move horizontally to the left and find the standardtopt value.(4)
Calculate the optimaltp time with the product between the standard topt value and α.The steps are schematically represented in Figures10–12.Figure 10
Step 1 for use of the proposed chart.Figure 11
Step 2 for use of the proposed chart.Figure 12
Step 3 for use of the proposed chart.Following Step 4, we find that the optimaltp value when performing the preventive activity is at tp=0.443⋅α=279.99h. Then, the PM policy is ready to be implemented. Notice that, for other values of β, an interpolation could be useful between the curves.


## 6. Discussion

It is interesting to analyze how the solution changes with altering one of the parameters of the problem. Sensitizing this to the shape parameterβ, for different values of p, a chart has been developed and shown in Figure 13.Figure 13
Sensitizing the standard optimal PM time (t) to different values of β, for different values of p.The unplotted area reveals the fact that it is not convenient to perform PM activities.In the case of sensitizing the values of the standard optimal PM timet for different values for the scale parameter α, considering the same β and p values, we find that there is a direct relation, and this has been shown in Section 4.From the results, it can be stated that if there are two kinds of equipment, Equipment A and Equipment B, both with the sameβ and p values, Equipment A with scale parameter αA, and Equipment B with scale parameter αB, and the optimal tp for Equipment A is tp,A, then the optimal tp for Equipment B is tp,B=αA/αB⋅tp,B.


## 7. Conclusion

From the present research, it is possible to conclude that, under PRP modeling and minimal long-run cost of maintenance in items, there is a boundary zone of possible combinations ofp and β parameters in which it is convenient to perform this kind of activity.Moreover, it has been shown from the solution for the optimization problem that a general solution could be obtainable by a chart that contains enough information to get the optimal PM time for those combinations ofp and β parameters, for which it is convenient to perform PM under PRP modeling. The solution to the problem is provided in the appendix.This solution has enormous potential because the tasks of the person in charge related to the determination of the policies of preventive maintenance or replacement are reduced drastically. Now, this person does not need the analytical skills and deep understanding of the mathematical problem that could cause considerable errors. This method indeed dramatically reduces the possibility of making calculation mistakes.Nowadays, increasing competitiveness is a requirement in the industry; this research is a contribution in terms of expediting processes and accelerated insertion of good policy implementation as the optimal preventive maintenance, and also in terms of education of reliability.Future research based on the elaboration of this type of instrument can be extended in general to other types of maintenance currently addressed in the literature, such as the consideration of imperfect maintenance, and, in turn, this tool can be applied and adapted for different models of failure processes that consider these characteristics, for example, NHPP and GRP models. Regarding the latter, the challenge will be to correctly manipulate and reduce the number of parameters that models include. Another area of future research addresses the development and implementation of generalized chart-based tools for maintenance models that consider a delay time stage, including in this way monitoring information on the equipment that conforms to the system to be analyzed.



---

*Source: 3078041-2020-09-24.xml*