Discrete model for a multi-objective maintenance optimization problem of safety systems

Citation: Briš, R.; Tran, N.T.T. Discrete Model for a Multi-Objective Maintenance Optimization Problem of Safety Systems. Mathematics 2023, 11, 320. https://doi.org/10.3390/ math11020320 Academic Editors: Michael Todinov and Ding Faxing Received: 4 November 2022 Revised: 18 December 2022 Accepted: 4 January 2023 Published: 7 January 2023 Copyright: © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). mathematics Article Discrete Model for a Multi-Objective Maintenance Optimization Problem of Safety Systems Radim Briš * and Nuong Thi Thuy Tran Department of Applied Mathematics, Faculty of Electrical Engineering and Computer Science, VSB—Technical University of Ostrava, 708 00 Ostrava-Poruba, Czech Republic *Correspondence: [email protected] Abstract: The aim of this article was to solve a multi-objective maintenance optimization problem by minimizing both unavailability and cost through the use of an optimal maintenance strategy. The problem took into account three different system designs upon which the objective functions are dependent, and the time to start preventive maintenance (PM) was used as a decision variable. This variable was optimized for all system components using a discrete maintenance model that allows for the specification of several discrete values of the decision variable in advance to find the optimal one. The optimization problem was solved using innovative computing methodology and newly updated software in MATLAB, which was used to quantify the unavailability of a complex system represented through a directed acyclic graph. A cost model was also developed to compute the cost of different maintenance configurations, and the optimal configuration was found. The results for a selected real system (a real fluid injection system adopted from references) showed that unavailability was less sensitive to variations in maintenance configurations, while cost variations were more noticeable in relation to different maintenance configurations. Applying PM, the increasing value of the decision variable increased cost because it led to more frequent corrective maintenance (CM) actions, and recovery times due to CM were more expensive than recovery times due to PM. Keywords: multi-objective optimization; unavailability; cost; maintenance; acyclic graph; alternating renewal process MSC: 60K10; 90B25 1. Introduction A real complex system can either be in a functioning or failed state. The probability of a functioning state under specified conditions over an intended period of time is usually defined [ 1 ] as system reliability R(t). System maintainability is defined to be the probability the system can be restored to a functional state within a specified period of time known as downtime. To model reliability and maintainability, we define two relevant random variables: time to system failure and time to repair. Maintainability and reliability are the two most important factors that must be considered while designing the system. They directly affect the availability A(t) of the system, which is defined as the probability that the system will continue to operate satisfactorily at any given point of time when it is being used under the specified conditions. The term reliability is often associated with systems that cannot be repaired; availability, however, is a term associated with repairable systems because it encompasses the full failure–recovery cycle over the mission time of the system. Thus, repairable systems operate during a time to failure (random variable) until a failure occurs. To recover the system’s operating state, time to repair (random variable) is needed. Both random variables are modelled by appropriate continuous probability distributions. Various maintenance strategies have been subjected to intense study and research in order to improve the reliability, availability, and usability of relevant industrial systems. Mathematics 2023,11, 320. https://doi.org/10.3390/math11020320 https://www.mdpi.com/journal/mathematics Mathematics 2023,11, 320 2 of 18 Unexpected failures can endanger human lives, cause unplanned production outages, etc., and this is why systems must be protected against them. A relevant tool to significantly improve system reliability is maintenance, both preventive and corrective. The authors in [ 2 ] state that maintenance is no longer a necessary evil and that production companies should invest in maintenance to maximize profit. To find the most suitable maintenance strategy for every component or subsystem, one must also consider the economic impact of the strategy employed to ensure profitable production. Therefore, different strategies must be evaluated in terms of their performance [2,3]. One of the common maintenance policies applied is a CM strategy, quite often denoted as a repair policy. It is launched at the time of failure when the system is already broken. Thereafter, the system is repaired into a functioning state through applying a CM action. The aim of a PM policy is to prevent the system from undesired breakdowns. PM is usually carried out while the system is still operational—it reduces the ageing processes, correlating to a decreased probability of system failure. Maintenance modeling is an emerging scientific discipline with rapid developments, and we do not want to present here a general overview of references on maintenance as this article is oriented towards the optimization of PM policy in context of system design. Several papers oriented towards both aspects are discussed below. Many techniques can be utilized to either maximize the system availability or minimize the unavailability. Unavailability is the complement of system availability to one. The present article is focused on two of them. First, this is done by modifying the design process for the system to ensure enough redundancies are present to reduce the system’s unavailability. For a system with a series-parallel configuration, redundancy can be thought of as a modification to the configuration to increase the number of parallel paths [ 4 ]. Second, an overall decrease in system unavailability is possible through PM [ 5 ]. The unavailability of a system due to its failure can occur at any time, requiring a significant endeavor to revert it back to the operating state. Contrarily, a planned shutdown to perform a PM task can represent a controlled situation with materials, spares, and human teams available, resulting in a reduced period of unavailability. It is clear now that the unavailability of repairable systems can be improved in two ways, either by modifying the design or applying PM strategies. Only a few research articles are devoted to looking at the simultaneous optimization of both ways from a multi-objective perspective. In [ 6 ], a methodology for integrated safety system design and maintenance optimization based on a bi-level evolutionary process was demonstrated. The authors tried to find both the optimum maintenance strategy and the optimum system design by applying genetic algorithms (GAs) as the optimization method and cost and unavailability as objective functions. The simultaneous optimization of design and maintenance during the life cycle was also presented in [ 7 ]. Optimization was performed using GAs, and objective functions were used for system reliability, redundancy, and life-cycle cost. A new approach to parallel optimization of maintenance and design of complex systems using reliability and cost as objective functions was demonstrated in [8]. The authors in [ 9 ] implemented optimization of both system design and PM strategy by coupling a multi-objective evolutionary algorithm and discrete simulation. System design can be optimized for reliability using redundant components, whereas PM strategy optimizes the PM times of each system component. The system availability and operation cost were the objective functions that were maximized and minimized, respectively. They applied a simulation approach in which each solution generated by the multi-objective evolutionary algorithm was evaluated through the use of a discrete simulation. This method explains the evolution of the system as it varies depending on operation and recovery times. This technique enables the analysis of complex real systems. Several configurations of the multi-objective evolutionary algorithm non-dominated sorting genetic algorithm II (NSGAII) were explored in [ 10 ]. The authors in [ 11 ] realized an exhaustive encoding comparative study, wherein some binary encoding alternatives were investigated. The authors were able to determine the optimal time to start a PM activity. All these approaches were based on a Mathematics 2023,11, 320 3 of 18 simulation technique that suffers from uncertainty, although the evolutionary processes can be improved in different ways, as for example, the effect of several chromosome lengths. The system unavailability and operation cost are the objective functions as well in this article, but their optimization was solved by an alternative method. The effective analytical method based on modelling was developed, proceeding from our previous findings in renewal theory and alternating renewal processes [ 12 ]. The theorem, called a recurrent linear integral equation, was modified to implement the new decision variable, i.e., time to start a PM activity, in the context of problem formulation. In addition, a new cost model for the system configuration corresponding to the desired PM strategy is defined in this paper. The innovative theorem, cost model, and optimization algorithm were numerically modeled using the high-performance programming language MATLAB. This article describes a new method to find the optimal PM policy to solve the designed optimization problem. The time to start PM was used as a decision variable in the optimization problem. This variable, which determines different maintenance modes of a system component, was optimally selected from a set of possible realistic maintenance modes. Optimization was performed for all system components. Thus, the discrete maintenance model was considered, wherein each component can function in one of several maintenance modes. The fixed value of the decision variable determined one maintenance mode of the component that predetermines both the evolution of unavailability and cost. Different maintenance modes of system components resulted in different system configurations, with each having a specific unavailability course as well as cost. The optimization process often demands plenty of computation time because a complex system can have several maintenance configurations. The discrete maintenance optimization is demonstrated on a real system selected from practice: the fluid injection system was adopted from the reference authorized by Cacereño et al. [11]. Literature on Maintenance Optimization Approaches There are different methods to solve the multi-objective maintenance optimization problem. A recent thorough classification can be found in [ 13 ]. In our article, we try to optimize the parameters of an a priori selected maintenance strategy. This problem can be solved by different approaches: • Our approach can be classified into mathematical approaches, wherein the optimization problem is formulated by means of mathematical equations, which are then solved by means of differential calculus to identify the optimal parameters of the maintenance strategy. In [ 14 ], a mathematical approach was used for optimizing maintenance profitability. Mathematical approaches were used in other references: for example, a single unit was optimized in [ 15 ] to optimize the scheduled maintenance strategy and the inventory management, optimal maintenance in the context of uncertainty is solved in [ 16 , 17 ], etc. Mathematical approaches can be used for the systems for which the optimization problem can be solved analytically or numerically. • Mixed integer programming is the field of optimization that addresses optimization problems with continuous and integer variables in the objective or in the constraints. Linear or nonlinear problems can be solved by means of the method in [ 18 ]. If the method is used for maintenance optimization, the possible maintenance optima are represented by integer variables. Examples of application of the method are the following: in [ 19 ], the authors optimized the maintenance schedule of a wind farm, in [ 20 ], a power distribution system was optimized; etc. In [ 21 ], the authors developed a mixed integer programming model for integration production and scheduled maintenance planning that considers the system’s manufacturing capacity and its reliability. The use of the method for maintenance optimization is mainly limited to simple systems because the computation time rapidly increases with the intricacy of systems [22]. • Dynamic programming is a method for solving multi-stage decision problems. The basic idea of the method is that complex problems are decomposed into simpler subproblems to be solved recursively at each time step. Examples of using the method Mathematics 2023,11, 320 4 of 18 are the following: the optimal maintenance strategy of road networks under budget constraints was solved in [ 23 ], the optimization of the maintenance check schedules in the aeronautical industry was determined in [ 24 ], the optimal maintenance strategy for power cables was solved in [ 25 ], etc. The main problems with dynamic programming are the curse of dimensionality and the need to explicitly define the transition probabilities among all the possible system states, which makes it inapplicable for complex systems [26]. • Metaheuristic search algorithms are computational processes where the solution of an optimization problem is found approximately by iteratively improving the candidate solutions [ 27 ]. For example, GAs are based on the principles of genetics and natural selection. GAs have been used in many situations to solve the maintenance optimization problems: the scheduled maintenance strategy of a wind farm was optimized in [28], the scheduled maintenance strategy of a multi-unit system was optimized in [ 29 ], the PM plan was optimized by means of multi-objective GAs in [ 30 ], etc. Moreover the maintenance optimization problems can be solved by other metaheuristic search algorithms: the particle swarm optimization algorithm was used to optimize the predictive maintenance interval of a manufacturing system [ 31 ], the harmony search algorithm was applied to find the best maintenance strategy for bridge infrastructures [ 32 ], simulated annealing was used to find the optimal scheduled maintenance plan of bridge networks [ 33 ], ant colony optimization was applied to optimize the maintenance scheduling of multi-unit systems [ 34 ], etc. These algorithms are easy to understand and easily adaptable to different optimization problems. A drawback of these algorithms is that they are slow to converge and do not guarantee convergence towards the global optimum. 2. Formulation of a Multi-Objective Optimization Problem Any multi-objective optimization problem works on the assumption that in general m-objective functions f 1 ( x ), f 2 ( x ), . . . f m ( x ), each one varying in a given range has to be optimized, i.e., either maximized or minimized, constrained by several restrictions imposed on the decision variables that are mostly related to system components. The optimization problem in this article is formulated using the two following objective functions: f 1 ( x ), which represents the cost function CS , and f 2 ( x ), which represents the stabilized unavailability function US . Thus, we searched for an optimal vector of decision variables x , here considered as the vector of times to start PM all of kcomponents (TP1, . . . , TPk) that minimize both f1(x) and f2(x): min [f1(x)∩f2(x)] (1) In our notation, we searched for the optimal decision vector (TP1, . . . , TPk) , minimizing both objective functions: min [CS(TP1, . . . , TPk)∩US(TP1, . . . , TPk)] (2) CS(TP1, . . . , TPk). . . total cost of maintenance of a system configuration; US(TP1, . . . , TPk). . . stabilized (asymptotic) system unavailability at the end of a mission time TM; (TP1, . . . , TPk). . . decision variable vector; k . . . number of system components, each having the decision variable TPi (time to start PM of i-th component), which is optimized in this article. From the point of view of computational feasibility, the optimization process was realized under following constraints: each decision variable has a prescribed domain containing three possible values, namely, TPmin , TPmax , and the middle point TPmid between minimal and maximal possible values. TPie[TPi,min;TPi,mid;TPi,max]i=1, . . . , k(3) Mathematics 2023,11, 320 5 of 18 Different values of TPi constitute different maintenance modes of the component, and for each maintenance mode, one can first compute the time evolution of unavailability function U(t) using the methodology described in Section 4.2, followed by the cost Ci,TM described in Section 4.3. Evolutions of U(t) of individual modes are then aggregated by means of the methodology based on directed acyclic graph (AG) described in Section 4.1 in order to compute the unavailability evolution of one system configuration. One system configuration is characterized by an asymptotic value US and total cost CS . All system configurations are finally ordered with respect to Equation (2). This is the main idea of the discrete maintenance model further introduced in Section 3. In most cases, both CS and US are usually complex linear or non-linear functions of the decision variable vector (TP1, . . . , TPk) , representing variables for which optimal values have to be found. Sometimes, the optimization problem is complicated by other constraints that apply for a given scope (e.g., a limitation for system performance). Such matters must be solved by finding the global optimum not violating any constraint. Thus, optimization, objective functions, and constraints cannot be managed independently because they influence each other. Solution of such optimization issues require advanced numerical algorithms [ 35 – 37 ]. The optimization problem with one objective function and restrictions including a specified limitation of maximal permissible value of the unavailability function was solved in our previous research work, for example, in [ 38 ]. In this article, we solved the optimization problem minimizing two objective functions (see Equation (2)) and respecting krestrictions given by Equation (3). Another tool to solve similar optimization problems results from GAs that have been frequently used in previous research, for example in [ 39 ], where GAs were used to optimize surveillance testing and maintenance, or in [ 11 ], which provided a real application example to demonstrate the innovative methodology presented in this article. As mentioned above, concerning the optimization outcome and decision variables, the discrete maintenance model can be classified as a process that finds optimized parameter values defining a single maintenance strategy selected a priori, e.g., in this paper, the time to start a PM activity must be found optimally, in the context of the problem formulation. To solve the multi-objective optimization problem, we selected a mathematical approach for the following reasons: • The search algorithm of the optimization problem requires good conditions for computing unavailability. We have the long-term experience to generate effective numerical algorithms for quantifying instantaneous unavailability of different types of maintained components, as well as complex systems. For example, in [ 12 ], we created and numerically elaborated the recurrent linear integral equation proceeding from alternating renewal processes. In this paper, the algorithms were further developed, i.e., properly modified and adopted to solve the formulated optimization problem. • We have long-term experience with the high-performance programming language MATLAB that was effectively used for the development of all numerical algorithms in the paper that are absolutely necessary to compute the optimization problem. • In the future, we intend to continue the research work in close collaboration with power industry experts, who are oriented towards the optimization of complex distribution networks that are hard to solve by the other above-mentioned alternative approaches. 3. Discrete Maintenance Model In [ 40 ], we introduced the discrete maintenance model for complex real systems with non-identical components. The model investigates systems with repairable components and latent failures to find optimal maintenance strategies to minimize cost and maintain unavailability under restriction. Latent failures are indicated using adaptable periods of inspections as a decision variable. A real complex system consists of a multitude of components that can be maintained in different ways, both by corrective and preventive interventions. Any component can Mathematics 2023,11, 320 6 of 18 function in different maintenance modes. One discrete maintenance mode of the ith component is determined by a selected value of the decision variable (in our application x i =TP i ) that immediately affects the maintenance cost of the mode. On the condition that a system consists of kcomponents, wherein each component can have five maintenance modes, in total, we have to explore 5 k maintenance configurations of the system. Each configuration is characterized by a typical unavailability value (mostly maximal or asymptotical) and total cost CS, which is commonly computed as a sum of the costs of all component modes constituting the configuration. It is necessary to find the optimal system configuration that meets requirements (2). We call this maintenance model the discrete maintenance model in this article. 4. Methodology for Computing Unavailability and Cost of a System Configuration 4.1. Graph Structure as a System Representation Any system structure can be represented with the help of a directed acyclic graph. Figure 1demonstrates the AGs of a real system from practice—a fluid injection system, which is later analyzed and optimized in detail in the next section. AGs are frequently used as systematic schemes to quantify the unavailability of complex systems [ 41 ] because they enable a reflective description of a system’s functionality. Obviously, AGs contain nodes and edges, wherein the exceptional node is the TOP node that describes the functionality of all systems depending on the functionality of its inferior subsystems and components that constitute internal and terminal nodes. Nodes are interconnected by edges, and AGs are acyclic, meaning that feedback loops are inadmissible. Terminal nodes—for example, V1 or P3—are denoted by blue squares, and they represent system components. Failure time as well as repair time or time necessary for PM of a system component are described by a suitable probability distribution. Applying these distributions, the time evolution of unavailability for each component can be computed using advanced renewal theory [ 40 ]. Mathematics 2023, 11, 320 7 of 18 Figure 1. Directed acyclic graph of a system from practice. 4.2. Model for Unavailability Exploration of a Terminal Node with Both CM and PM To find the unavailability function of TOP node U(t), it is necessary to find a model and algorithm for the unavailability quantification of terminal nodes that undergo both PM and CM. At first, the model with CM will be introduced, which will be further generalized to enable the implementation of both PM and CM. Applying CM, we have to consider two mutually cooperating random variables: the lifetime X, described by either distribution function F(t) or probability density function (pdf) f(t), and a random time necessary for completing the repair, actually repair, or recovery time Y, described by either distribution function G(t) or pdf g(t). Resulting from renewal theory and alternating renewal processes, availability A(t) can be computed as follows [12]: 𝐴(𝑡) =1−𝐹(𝑡)+∫ℎ(𝑥)[1−𝐹(𝑡−𝑥)]dx 𝑡 0=𝑅(𝑡)+∫ℎ(𝑥)𝑅(𝑡−𝑥)dx, 𝑡 0 (4) where R(t) = 1 − F(t) is the reliability function and h(x) is the renewal density of the corresponding alternating renewal process. The unavailability U(t) is given as follows:   −−−=−= t dxxtFxhtFtAtU 0 )(1)()()(1)( (5) To compute the unavailability from the formula (5), h(x) must be known, which could be a problem in practice, bringing about a large number of numerical complications because the renewal density is numerically represented as an infinite sum of probability densities—here, each being computed as a convolution. Fortunately, formula (5) can be superseded by the equivalent formula (6), which brings the following theorem, called the recurrent linear integral equation, which was first mentioned and proven in [12]. Theorem 1. The unavailability U(t) in formula (5) is equivalent to the unavailability U(t) in the following formula (6) Figure 1. Directed acyclic graph of a system from practice. Internal nodes are denoted by blue triangles, e.g., u1 or u2, respresenting subsystems of the fluid system. Subsystems as well as components (i.e., internal respectively terminal nodes) at a given time are either correctly working or in a failed state (under restoration). A node is in a functioning state when the number of child nodes is greater or equal to the Mathematics 2023,11, 320 7 of 18 number of nodes inside the triangle, otherwise it is in a failed state. For example, the node u1 is correctly working if the number of correctly working directly inferior nodes is either 1 or 2. The knowledge of the unavailability functions of terminal nodes (components) can be used to quantify the unavailability function of internal nodes (subsystems), with both serving as inputs for computing the unavailability function of the TOP node U(t), which represents all systems. The unavailability function U(t) demonstrates the time-dependent probability that the system is unavailable at time tdue to a failure or due to a still ongoing repair process. 4.2. Model for Unavailability Exploration of a Terminal Node with Both CM and PM To find the unavailability function of TOP node U(t), it is necessary to find a model and algorithm for the unavailability quantification of terminal nodes that undergo both PM and CM. At first, the model with CM will be introduced, which will be further generalized to enable the implementation of both PM and CM. Applying CM, we have to consider two mutually cooperating random variables: the lifetime X, described by either distribution function F(t) or probability density function (pdf) f(t), and a random time necessary for completing the repair, actually repair, or recovery time Y, described by either distribution function G(t) or pdf g(t). Resulting from renewal theory and alternating renewal processes, availability A(t) can be computed as follows [ 12 ]: A(t) = 1−F(t) + Zt 0h(x)[1−F(t−x)]dx =R(t) + Zt 0h(x)R(t−x)dx, (4) where R(t) = 1 − F(t) is the reliability function and h(x) is the renewal density of the corresponding alternating renewal process. The unavailability U(t) is given as follows: U(t) = 1−A(t) = F(t)− t R 0 h(x)[1−F(t−x)]dx (5) To compute the unavailability from the Formula (5), h(x) must be known, which could be a problem in practice, bringing about a large number of numerical complications because the renewal density is numerically represented as an infinite sum of probability densities—here, each being computed as a convolution. Fortunately, Formula (5) can be superseded by the equivalent Formula (6), which brings the following theorem, called the recurrent linear integral equation, which was first mentioned and proven in [12]. Theorem 1. The unavailability U(t) in Formula (5) is equivalent to the unavailability U(t) in the following Formula (6) U(t) = Zt 0f(x)[1−G(t−x)]dx +Zt 0(f∗g)(x)U(t−x)dx (6) where ∗means convolution. Following this, we consider the PM strategies to maintain the operating status of the node. These operations performed before failure are usually less complex than CM activities that must be undertaken in the case of failure. Obviously, each PM activity starting at time TP that serves as a decision variable in our discrete maintenance model requires some recovery time, which we define to be a new random variable Z. Thus, the operational time of a component is interrupted either by the time to failure X(lifetime) or by the time to start a PM activity TP, whichever occurs first. A recovery time can be realized either by the recovery time Ydue to CM or recovery time Zdue to PM. In other words, PM activities are scheduled shutdowns, and recovery times are usually shorter and less expensive than Mathematics 2023,11, 320 8 of 18 the same due to CM. Examples of PM include ensuring the availability of spare parts and the training of human personnel. PM activities should be conducted optimally before the failure but as close as possible to it. From an optimization point of view, it is necessary to minimize the system unavailability as well as cost due to recovery times. To compute component unavailability U(t) respecting both CM and implemented PM activities, our Formula (6) can be employed with the modification that the random variable Xwill be substituted by the random variable V= min(X,TP), which represents interruption of the operation time. The distribution function FVof Vcan be easily found: FV(t) = P(min(X,TP)<t) = P(X<t) + P(TP <t)−P(X<t). P(TP <t) = 1−P(min(X,TP)≥t)= 1−P(X≥t). P(TP ≥t)(7) Consequently, the following formulas hold: FV(t) = F(t)for t<TP FV(t) = 1 for t≥TP (8) Thus, in the latter case, theorem (6) is modified as follows: U(t) = TP Z 0 f(x)·[1−G(t−x)]dx+ TP Z 0 (f∗g)(x)·U(t−x)dx + t−TP Z TP w(x)U(t−TP −x)dx (9) where the last integral of (9) mathematically represents the remaining contribution to unavailability function U(t) for t ≥ TP, provided that at time TP the recovery time started due to PM. w(x) is the pdf of recovery time Z. The expected value µVof Vcan be found according to the following formula: E V =µV=ZTP 0(1−FV(t))dt (10) Computing possibilities of the method for unavailability quantification of complex multi-component and highly reliable systems were successfully demonstrated in a comparison study in [42], as well as in [12]. 4.3. Cost Model of a System Configuration The cost model of a system configuration can be obtained by adding up all contributions resulting from both CM and PM replacement interventions of a mode over all of the system components. Components can operate in different maintenance modes; the cost of one maintenance mode consists of two main contributions generated by CM on the one hand and PM on the other. The cost of CM further depends on the mean number of all recovery times due to both CM and PM during mission time T M and CM parameters. The cost of PM depends on the decision variable TP and PM’s parameters. In practical situations, the cost contributions result from a year database to gain an average yearly cost for system configurations in a monitored period. In the remainder of this article, the cost computed in non-identified cost units on the basis of the summation principle is provided. To obtain the cost of one system configuration, we simply add up the costs of all maintenance modes of all system components. The cost of one maintenance mode of i-th component Ci,TMcan be computed as follows: Ci,TM=ni,R.F(TPi).Ci,R+ni,R.R(TPi).Ci,PM (11) where ni,R=TM MTTIi+MRTi. . . the mean number of recovery actions of the i-th component per mission TM; MTTIi=µV. . . is the mean time to intervention caused by either CM or PM; MRTi. . . mean recovery time of the i-th component, which is due to either PM or CM. Mathematics 2023,11, 320 9 of 18 MRTi=F(TPi).MTRCMi+R(TPi).MTRPMi(12) MTRCMi. . . mean recovery time due to CM; MTRPMi. . . mean recovery time due to PM; TPi. . . decision variable of the i-th component determining PM strategy; Ci,R. . . CM cost = cost of one CM intervention of the i-th component in cost units; Ci,PM . . . PM cost = cost of one PM intervention of the i-th component in cost units. The total cost of one system configuration CS is given by summing up these contributions described by Formula (11) over all of the system components k: CS=∑k i=1Ci,TM(13) 5. Results with the Real Complex System—Fluid Injection System and Discussion The discrete maintenance optimization was demonstrated on the industrial fluid injection system adopted by authors in [11], wherein the authors created a massive simulation to study the system using GA to achieve both maximum availability and minimum cost, paying attention to the possible impact on solutions as a result of different encodings, chromosome lengths, and parameter configurations. This article differs in computing methodology—it is based on probabilistic modelling resulting from alternating renewal processes. It was necessary to optimize both the design and the PM policy, i.e., to find the optimal decision variable vector x = (TP 1 , . . . ,TP k ) minimizing both unavailability and cost. The system consisted of valves (V) and pumps (P) and is depicted in its full version in Figure 2. Two other partial versions (designs) were taken into account: Version 1 intended as the simplest version without parallel ordered components P2 and V4, and Version 2, which was intended as the full version without pump P2. The functionality of the system (full version) is described and analyzed within the framework of Section 4.1, utilizing the AG depicted in Figure 1. Unavailability of the system can be decreased by increasing investment in the PM on the one hand, but it signifies the growth of cost on the other hand. Mathematics 2023, 11, 320 9 of 18 system configurations in a monitored period. In the remainder of this article, the cost computed in non-identified cost units on the basis of the summation principle is provided. To obtain the cost of one system configuration, we simply add up the costs of all maintenance modes of all system components. The cost of one maintenance mode of i-th component 𝐶𝑖,𝑇𝑀 can be computed as follows: 𝐶𝑖,𝑇𝑀=𝑛𝑖,𝑅.𝐹(𝑇𝑃𝑖).𝐶𝑖,𝑅 +𝑛𝑖,𝑅.𝑅(𝑇𝑃𝑖).𝐶𝑖,𝑃𝑀 (11) where 𝑛𝑖,𝑅 =𝑇𝑀 𝑀𝑇𝑇𝐼𝑖+𝑀𝑅𝑇𝑖 … the mean number of recovery actions of the i-th component per mission TM; 𝑀𝑇𝑇𝐼𝑖 = μV …is the mean time to intervention caused by either CM or PM; 𝑀𝑅𝑇𝑖 … mean recovery time of the i-th component, which is due to either PM or CM. 𝑀𝑅𝑇𝑖=𝐹(𝑇𝑃𝑖).𝑀𝑇𝑅𝐶𝑀𝑖+𝑅(𝑇𝑃𝑖).𝑀𝑇𝑅𝑃𝑀𝑖 (12) 𝑀𝑇𝑅𝐶𝑀𝑖 … mean recovery time due to CM; 𝑀𝑇𝑅𝑃𝑀𝑖 … mean recovery time due to PM; TPi …decision variable of the i-th component determining PM strategy; 𝐶𝑖,𝑅 … CM cost = cost of one CM intervention of the i-th component in cost units; 𝐶𝑖,𝑃𝑀 …PM cost = cost of one PM intervention of the i-th component in cost units. The total cost of one system configuration 𝐶𝑆 is given by summing up these contributions described by formula (11) over all of the system components k: 𝐶𝑆=∑ 𝐶𝑖,𝑇𝑀 𝑘 𝑖=1 (13) 5. Results with the Real Complex System—Fluid Injection System and Discussion The discrete maintenance optimization was demonstrated on the industrial fluid injection system adopted by authors in [11], wherein the authors created a massive simulation to study the system using GA to achieve both maximum availability and minimum cost, paying attention to the possible impact on solutions as a result of different encodings, chromosome lengths, and parameter configurations. This article differs in computing methodology—it is based on probabilistic modelling resulting from alternating renewal processes. It was necessary to optimize both the design and the PM policy, i.e., to find the optimal decision variable vector x = (TP1, …, TPk) minimizing both unavailability and cost. The system consisted of valves (V) and pumps (P) and is depicted in its full version in Figure 2. Two other partial versions (designs) were taken into account: Version 1 intended as the simplest version without parallel ordered components P2 and V4, and Version 2, which was intended as the full version without pump P2. The functionality of the system (full version) is described and analyzed within the framework of Section 4.1, utilizing the AG depicted in Figure 1. Unavailability of the system can be decreased by increasing investment in the PM on the one hand, but it signifies the growth of cost on the other hand. Figure 2. Real fluid injection system adopted from [11]. Figure 2. Real fluid injection system adopted from [11]. Other assumptions: •only two component states are admissible: operational and faulty state; •the components are mutually independent; •if failure of a component comes, the repair starts immediately; • if a repair is completed, the component’s state is equivalent to that of a new component. Table 1brings about the data used in the analysis. Additional information on the input parameters is found in the Notations. The data were adopted from the special source for reliability analysis OREDA (offshore reliability data handbook [ 43 ]), further from expert appraisal (originating from the experience of the Machinery and Reliability Institute (MRI), Alabama, USA) and corresponding reliability mathematics. The OREDA has traditionally been focused on the reliability and availability of production systems and equipment. As we declared above, the discrete maintenance optimization process was carried out, resulting from the minimization of both unavailability and cost due to maintaining system stages including both PM and CM. To realize this, one must Mathematics 2023,11, 320 16 of 18 Notations Variables Xtime to failure (the lifetime) Yrepair (recovery) time after a failure occurs Ztime to perform a PM TP decision variable—time to start a PM process V= min(X,TP) time to interruption of the operation time caused by either failure or PM (TP1, . . . , TPk)vector of decision variables of kcomponents, days TPidecision variable of i-th component determining its PM strategy, days TPi,min;TPi,max minimal and maximal possible values of TPi TPi,mid TPi,min +TPi,max/2 CS(TP1, . . . , TPk)total cost of maintenance of a system configuration, cost units US(TP1, . . . , TPk)stabilized (asymptotic) system unavailability at the end of TM ni,Rmean number of recovery actions of i-th component per mission time TM MTTIimean time to intervention of i-th component caused by either CM or PM, days MRTimean recovery time of i-th component that is due to either PM or CM, days MTRCMimean recovery time of i-th component due to CM, days MTRPMimean recovery time i-th component due to PM, days Indices F(t) distribution function of a random variable X R(t)=1−F(t) reliability function of a random variable X f(t) probability density function (pdf) of a random variable X G(t) distribution function of a random variable Y g(t) probability density function (pdf) of a random variable Y W(t) distribution function of a random variable Z w(t) probability density function (pdf) of a random variable Z FV(t) distribution function of a random variable V U(t) instantaneous time-dependent unavailability function A(t)=1−U(t) instantaneous availability function h(x) renewal density Parameters Ci,RCM cost = cost of one CM intervention of the i-th component in cost units Ci,PM PM cost = cost of one PM intervention of the i-th component in cost units TMmission time, days λthe failure rate: parameter of exponential distribution of the random variable X, h−1 [TRmin;TRmax] parameters of rectangular distribution of the random variable Y [TRPmin;TRPmax] parameters of rectangular distribution of the random variable Z References 1. Misra, K.B. Reliability Engineering: A Perspective. In Handbook of Performability Engineering; Misra, K.B., Ed.; Springer: London, UK, 2008; pp. 253–289. [CrossRef] 2. Ding, S.-H.; Kamaruddin, S. Maintenance policy optimization—Literature review and directions. Int. J. Adv. Manuf. Technol. 2015 , 76, 1263–1283. [CrossRef] 3. Lee, J.; Lapira, E.; Bagheri, B.; Kao, H.-A. Recent advances and trends in predictive manufacturing systems in big data environment. Manuf. Lett. 2013,1, 38–41. [CrossRef] 4. Kuo, W.; Prasad, R.; Tillman, F.A.; Mwang, C.L. Optimal Reliability Design: Fundamentals and Applications, 1st ed.; Cambridge University Press: Cambridge, UK, 2001; p. 389. 5. Gao, Y.; Feng, Y.; Zhang, Z.; Tan, J. An optimal dynamic interval preventive maintenance scheduling for series systems. Reliab. Eng. Syst. Saf. 2015,142, 19–30. [CrossRef] 6. Galván, B.; Winter, G.; Greiner, D.; Salazar, D.; Méndez, M. New Evolutionary Methodologies for Integrated Safety System Design and Maintenance Optimization. In Computational Intelligence in Reliability Engineering: Evolutionary Techniques in Reliability Analysis and Optimization; Levitin, G., Ed.; Springer: Berlin/Heidelberg, Germany, 2007; pp. 151–190. [CrossRef] 7. Okasha, N.M.; Frangopol, D.M. Lifetime-oriented multi-objective optimization of structural maintenance considering system reliability, redundancy and life-cycle cost using GA. Struct. Saf. 2009,31, 460–474. [CrossRef] 8. Adjoul, O.; Benfriha, K.; El Zant, C.; Aoussat, A. Algorithmic Strategy for Simultaneous Optimization of Design and Maintenance of Multi-Component Industrial Systems. Reliab. Eng. Syst. Saf. 2021,208, 107364. [CrossRef] Mathematics 2023,11, 320 17 of 18 9. Cacereño, A.; Greiner, D.; Galván, B. Solving Multi-objective Optimal Design and Maintenance for Systems Based on Calendar Times Using NSGA-II. In Advances in Evolutionary and Deterministic Methods for Design, Optimization and Control in Engineering and Sciences; Gaspar-Cunha, A., Periaux, J., Giannakoglou, K.C., Gauger, N.R., Quagliarella, D., Greiner, D., Eds.; Springer Nature: Cham, Switzerland, 2021; pp. 245–259. [CrossRef] 10. Deb, K.; Pratap, A.; Agarwal, S.; Meyarivan, T. A fast and elitist multiobjective genetic algorithm: NSGA-II. IEEE Trans. Evol. Comput. 2002,6, 182–197. [CrossRef] 11. Cacereño, A.; Greiner, D.; Galván, B.J. Multi-Objective Optimum Design and Maintenance of Safety Systems: An In-Depth Comparison Study Including Encoding and Scheduling Aspects with NSGA-II. Mathematics 2021,9, 1751. [CrossRef] 12. Briš, R.; Byczanski, P. On innovative stochastic renewal process models for exact unavailability quantification of highly reliable systems. Proc. Inst. Mech. Eng. Part O J. Risk Reliab. 2017,231, 617–627. [CrossRef] 13. Pinciroli, L.; Baraldi, P.; Zio, E. Maintenance optimization in Industry 4.0. Reliab. Eng. Syst. Saf.. Manuscript Draft JRESS-D-2201170R1 in press. 14. Oke, S.A. An analytical model for the optimisation of maintenance profitability. Int. J. Prod. Perform. Manag. 2005 ,54, 113–136. [CrossRef] 15. Rezg, N.; Chelbi, A.; Xie, X. Modeling and optimizing a joint inventory control and preventive maintenance strategy for a randomly failing production unit: Analytical and simulation approaches. Int. J. Comput. Integr. Manuf. 2005 ,18, 225–235. [CrossRef] 16. de Jonge, B.; Klingenberg, W.; Teunter, R.; Tinga, T. Optimum maintenance strategy under uncertainty in the lifetime distribution. Reliab. Eng. Syst. Saf. 2015,133, 59–67. [CrossRef] 17. Compare, M.; Martini, F.; Zio, E. Genetic algorithms for condition-based maintenance optimization under uncertainty. Eur. J. Oper. Res. 2015,244, 611–623. [CrossRef] 18. Sahinidis, N.V. Mixed-integer nonlinear programming 2018. Optim. Eng. 2019,20, 301–306. [CrossRef] 19. Besnard, F.; Patrikssont, M.; Strombergt, A.-B.; Wojciechowskit, A.; Bertling, L. An optimization framework for opportunistic maintenance of offshore wind power system. In Proceedings of the 2009 IEEE Bucharest PowerTech, Bucharest, Romania, 28 June–2 July 2009; pp. 1–7. [CrossRef] 20. Dehghani, N.L.; Darestani, Y.M.; Shafieezadeh, A. Optimal Life-Cycle Resilience Enhancement of Aging Power Distribution Systems: A MINLP-Based Preventive Maintenance Planning. IEEE Access 2020,8, 22324–22334. [CrossRef] 21. Aghezzaf, E.-H.; Khatab, A.; Le Tam, P. Optimizing production and imperfect preventive maintenance planning’s integration in failure-prone manufacturing systems. Reliab. Eng. Syst. Saf. 2016,145, 190–198. [CrossRef] 22. Putz, D.; Schwabeneder, D.; Auer, H.; Fina, B. A comparison between mixed-integer linear programming and dynamic programming with state prediction as novelty for solving unit commitment. Int. J. Electr. Power Energy Syst. 2020 ,125, 106426. [CrossRef] 23. Ma, J.; Cheng, L.; Li, D. Road Maintenance Optimization Model Based on Dynamic Programming in Urban Traffic Network. J. Adv. Transp. 2018,2018, 4539324. [CrossRef] 24. Deng, Q.; Santos, B.F.; Curran, R. A practical dynamic programming based methodology for aircraft maintenance check scheduling optimization. Eur. J. Oper. Res. 2020,281, 256–273. [CrossRef] 25. Sachan, S.; Zhou, C. Probabilistic dynamic programming algorithm: A solution for optimal maintenance policy for power cables. Life Cycle Reliab. Saf. Eng. 2019,8, 117–127. [CrossRef] 26. Sutton, R.S.; Barto, A. Reinforcement Learning: An Introduction, 2nd ed.; MIT Press: Cambridge, MA, USA, 2018. 27. Nesmachnow, S. An overview of metaheuristics: Accurate and efficient methods for optimisation. Int. J. Metaheuristics 2014 , 3, 320. [CrossRef] 28. Carlos, S.; Sánchez, A.; Martorell, S.; Marton, I. Onshore wind farms maintenance optimization using a stochastic model. Math. Comput. Model. 2013,57, 1884–1890. [CrossRef] 29. Dao, C.D.; Zuo, M.J. Selective maintenance of multi-state systems with structural dependence. Reliab. Eng. Syst. Saf. 2017 ,159, 184–195. [CrossRef] 30. Yulan, J.; Zuhua, J.; Wenrui, H. Multi-objective integrated optimization research on preventive maintenance planning and production scheduling for a single machine. Int. J. Adv. Manuf. Technol. 2008,39, 954–964. [CrossRef] 31. Loganathan, M.K.; Gandhi, O.P. Maintenance cost minimization of manufacturing systems using PSO under reliability constraint. Int. J. Syst. Assur. Eng. Manag. 2016,7, 47–61. [CrossRef] 32. García-Segura, T.; Yepes, V.; Frangopol, D.M.; Yang, D.Y. Lifetime reliability-based optimization of post-tensioned box-girder bridges. Eng. Struct. 2017,145, 381–391. [CrossRef] 33. Mao, X.; Jiang, X.; Yuan, C.; Zhou, J. Modeling the Optimal Maintenance Scheduling Strategy for Bridge Networks. Appl. Sci. 2020,10, 498. [CrossRef] 34. Tran, L.V.; Huynh, B.H.; Akhtar, H. Ant Colony Optimization Algorithm for Maintenance, Repair and Overhaul Scheduling Optimization in the Context of Industrie 4.0. Appl. Sci. 2019,9, 4815. [CrossRef] 35. Fletcher, R. Practical Methods of Optimization; Wiley: New York, NY, USA, 1987. 36. Harunuzzaman, M.; Aldemir, T. Optimization of Standby Safety System Maintenance Schedules in Nuclear Power Plants. Nucl. Technol. 1996,113, 354–367. [CrossRef] Mathematics 2023,11, 320 18 of 18 37. Vaurio, J.K. Optimization of test and maintenance intervals based on risk and cost. Reliab. Eng. Syst. Saf. 1995 ,49, 23–36. [CrossRef] 38. Briš, R.; Tran, N.T.T. Optimization of maintenance policies for complex and highly reliable multi-unit systems. In Safety and Reliability—Theory and Applications;ˇ Cepin, M., Briš, R., Eds.; Taylor & Francis Group: London, UK, 2017; pp. 403–411. 39. Munõz, A.; Martorell, S.; Serradell, V. Genetic algorithms in optimizing surveillance and maintenance of components. Reliab. Eng. Syst. Saf. 1997,57, 107–120. [CrossRef] 40. Briš, R.; Byczanski, P.; Goˇno, R.; Rusek, S. Discrete maintenance optimization of complex multi-component systems. Reliab. Eng. Syst. Saf. 2017,168, 80–89. [CrossRef] 41. Briš, R. Parallel simulation algorithm for maintenance optimization based on directed Acyclic Graph. Reliab. Eng. Syst. Saf. 2008 , 93, 874–884. [CrossRef] 42. Briš, R.; Tran, N.T.T. Newly enhanced computing algorithm to quantify unavailability of maintained multi-component systems. In Safety and Reliability—Safe Societies in a Changing World; Haugen, S., Barros, A., van Gulijk, C., Kongsvik, T., Vinnem, J.E., Eds.; CRC Press: Boca Raton, FL, USA, 2018; pp. 931–936. ISBN 9780815386827. 43. Offshore and Onshore Reliability Data. Offshore Reliability Data Handbook, 5th ed.; Stiftelsen for Industriell og Teknisk Forskning: Trondheim, Norway; Det Norske Veritas: Hovik, Norway, 2009; ISBN 9788214048308. Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.