Applications of artificial intelligence for ... - Semantic Scholar
ARTICLE IN PRESS
Engineering Applications of Artificial Intelligence 19 (2006) 113–126 www.elsevier.com/locate/engappai
Applications of artificial intelligence for optimization of compressor scheduling Hanh H. Nguyen, Christine W. Chan Faculty of Engineering, University of Regina, Regina, Saskatchewan, Canada S4S 0A2 Received 2 April 2005; received in revised form 13 June 2005; accepted 16 June 2005 Available online 10 October 2005
Abstract This paper presents a feasibility study of evolutionary scheduling for gas pipeline operations. The problem is complex because of several constraints that must be taken into consideration during the optimization process. The objective of gas pipeline operations is to transfer sufficient gas from gas stations to consumers so as to satisfy customer demand with minimum costs. The scheduling involves selection of a set of compressors to operate during a shift. The scheduling decision has to be made so as to satisfy the dual objectives of minimizing the sum of fuel cost, start-up cost, the cost of gas wasted due to oversupply, and satisfying minimal operative and inoperative time of the compressors. The problem was decomposed into the two subproblems of gas load forecast and selection of compressors. Neural networks were used for forecasting the load; and genetic algorithms were used to search for a near optimal combination of compressors. The study was conducted on a subsystem of the pipeline network located in southeastern Saskatchewan, Canada. The results are compared with the solutions generated by an expert system and a fuzzy linear programming model. r 2005 Elsevier Ltd. All rights reserved. Keywords: Genetic algorithms; Neural networks; Compressor scheduling
1. Introduction Evolutionary scheduling for gas pipeline compressor is a complex problem because of constraints that must be taken into consideration during the optimization process. Such constraints include operation cost, maintenance cost, and customer satisfaction. This paper reports on a feasibility study on applying Artificial Intelligence technologies to natural gas pipeline operations. The focus of the study was on the St. Louis East gas pipeline distribution system, a subsystem of the pipeline network located in southeastern Saskatchewan. The system consists of two stations located at Melfort and St. Louis. The two stations use compressors to Corresponding author. Tel.: +1 306 585 5225; fax: +1 306 585 4855. E-mail addresses: [email protected] (H.H. Nguyen), [email protected] (C.W. Chan).
0952-1976/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.engappai.2005.06.008
generate sufficient pressure to transport gas to the surrounding consumption areas of St. Brieux, Nipawin and Hudson Bay. The dispatchers’ task is to adjust the compression level in order to generate the necessary pressure while not wasting energy. When the customer demand increases, a dispatcher adds compression to the pipeline system by turning on one or more compressors. On the other hand, the dispatcher turns off one or more compressors to reduce compression in the pipeline system when the customer demand decreases. This decision has a significant impact on effectiveness of the natural gas pipeline operation as well as on operating costs. The objective of this study is to use automation techniques to aid the dispatcher in optimizing natural gas pipeline operations in order to satisfy customer demand with minimal operating costs. A dispatcher needs to know in advance when a large volume of gas is required and what way is the optimal to ensure it is
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
available. Previous studies were conducted in which expert system and linear programming (LP) technologies were employed to solve the problem. Here, we present a genetic algorithm (GA) system that generates workable schedules of compressor operation based on the current situation. Since the problem can be decomposed into two sub-problems of demand forecast and compressor selection, neural networks are employed to tackle the former problem of demand forecast and GAs are used for addressing the problem of compressor optimization. This paper is organized as follows. Section 2 provides some background on the past solutions for the problem and the GA technique in general. Section 3 describes the first subtask, as well as development and results of demand prediction using neural networks. Section 4 discusses the second subtask, which involves analysis and implementation of GAs for optimization of compressor scheduling. Section 5 presents comparison of results generated by the three methods. Section 6 gives some conclusion and discusses some possible future work.
2. Background 2.1. Problem domain This project was part of a study on automation of natural gas pipeline operations conducted jointly with a gas pipeline company in Saskatchewan, Canada. The objective of the study was on optimization of compressor scheduling and the study focused on a section of the Saskatchewan gas pipeline around the St. Louis East compressor station in the province of Saskatchewan, Canada. A schematic of the St. Louis East system is shown in Fig. 1. The system consists of two compressor stations, Melfort and St. Louis. These compressor stations supply natural gas to two customer locations, Nipawin and Hudson Bay. In St. Louis, there are three compressor units; two of them are electrical compressor units and
Melfort compressor station
Nipawin
Hudson Bay Electricity compressor Gas compressor
Fig. 1. Schematic of the St. Louis East system.
Industrial customer Time
2
Pressure
1
Pressure
St. Louis compressor station
the other is a gas compressor unit. In Melfort, there are two gas compressors. An electrical compressor unit provides 250 hp and a gas unit provides 600 or 650 hp. The demand for natural gas from customers fluctuates depending on the season. In the winter, the demand for natural gas is usually higher than in the summer. In addition, the demand for natural gas also changes depending on the time of day. For example, in the morning, the demand is higher because the customer begins to use natural gas to warm the premise. In the afternoon, the demand is low since the facilities are already heated up in the morning. The customers can also be grouped into three types: industrial, dehydrator, and heat sensitive customers. Each type of customer reflects a different pattern of natural gas consumption as illustrated in Fig. 2. The industrial customer has the same rate of natural gas consumption any time of the day. The dehydrator customer demands two set amounts of natural gas at different periods of time. For example, between 8:00 and 10:00 a.m., the demand is 240 103 m3/day while between 10:00 am to 12:00 am, the demand is 200 103 m3/day. The demand of the heat sensitive customer fluctuates over time and is difficult to predict. In practice, the demand for natural gas can fluctuate from 200 103 m3/day to over 560 103 m3/day within 1 or 2 h. A demand of 200 103 m3/day can be handled from St. Louis with one compressor unit. When the demand exceeds 560 103 m3/day, all units at St. Louis and both units at Melfort are needed. To ensure customer satisfaction, the operator needs to know beforehand and be ready for the largest volume requirement of gas. Otherwise, the system pressures at Nipawin and Hudson Bay will be below the required minimum. The graphs showing supply and demand of gas in the pipelines are generated during operations by the simulator program of the gas pipeline transmission company. The display includes a compressor discharge pressure curve vs. customer station pressure curve as shown in Fig. 3. The two curves in Fig. 3 provide pressure information to dispatchers, and indicate the relationship between demand and supply. The gap between the two curves is called the comfort zone (CZ). If this CZ is wide, then customer satisfaction is
Dehydrator customer Time
3
Pressure
114
Heat sensitive customer Time
Fig. 2. Different natural gas consumption patterns.
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126 6500 6000
Pressure (kpa)
5500 5000 4500 4000 3500 3000
Compressor Discharge Pressure (supply) Customer End-point Delivery Pressure (demand)
2500 0
200
400
600
800
Time (hr)
Fig. 3. Comfort zone.
guaranteed but the cost is high; if it is too narrow, then the operation is cost effective but future customer demand may not be satisfied. Hence, the dispatcher needs to monitor the volume of gas supplied to the customer and maintain a suitable width of the CZ between the two pressure curves. Optimizing compressor scheduling and deciding which of the five compressors to turn on/off is a primary step for increasing pressure in the pipelines while maintaining the optimum operating pressure in the pipeline system. A novice dispatcher may make unnecessary start/stop decisions due to lack of experience. A minimum of 4 months is required to train new dispatchers to operate the system. However, it takes a number of years for them to gain proficiency in operating the system smoothly and cost effectively. 2.2. Compressor scheduling: expert system and linear programming approaches The problem of compressor schedule optimization has been previously tackled with the technologies of expert system (Uraikul et al., 2000) and fuzzy LP (Sun et al., 2000). The expert system tackled the three subtasks of (1) determining the state of the current line pack to decide if compression needs to be added or reduced, (2) determining the required break horsepower (BHP), and (3) selecting the combination of compressors from the two stations to meet the BHP requirement with minimal cost. Since the main emphasis of the expert system was to support the user in the first two tasks of determining the control action and BHP, a rather simplistic approach was adopted for the third task of compressor selection. The system only took into consideration the BHP requirement as the basis for scheduling compressors. The system recommended five combinations of compressors to be turned on in the following order: (1) Free flow: St. Louis E2, (2) 0oBHPo800: St. Louis E2 and N1, (3) 800pBHP o1200: St. Louis E2 and N1,
115
and Melfort N2, (4) 1200pBHP o1600: St. Louis E2, E1 and N1, Melfort N2, and (5) 1600pBHPo2200: St. Louis E2, E1, N1, and Melfort N2 and N3, where St. Louis E1 and E2 are electrical compressors, and St. Louis N1, Melfort N2 and N3 are natural gas compressors. This suggestion is made based on heuristic knowledge of senior pipeline operators. The expert system has been demonstrated using simulation data for various operational scenarios to a panel of expert operators, and they have confirmed that its recommendations were valid. However, no attempts have been made to quantify the performance in terms of component costs. While customer demand satisfaction was an implicit assumption in the experts’ formulation of the rules in the expert system, other component costs were not explicitly considered. As a result, neither the rules nor the validation results reveal how well the expert system handles the other costs. A closer examination on these rules reveals some interesting hidden heuristics. For example, the rules involve only either addition or substraction of compressors, not both at the same time. Given the fact that it is unlikely for BHP requirements to oscillate between two ranges of BHP in a short time period of 3 h or less, the chance of incurring wear and tear on the units is limited. In terms of cost, a disadvantage of the expert system approach is that it enforces safety by sacrificing energy efficiency. The expert system rules are constructed in a way that ensures the upper bound of every range of BHP is always satisfied. Therefore, the rule on compressor selection for that range is not efficient when the actual BHP requirement is in fact closer to the lower bound of the range. A weakness of the expert system is that its decisions are based solely on data from instaneous snap-shots of the pipeline and hence ignore some important information such as demand trend. On the other hand, the LP and GA approaches take into consideration data from the entire duration of a shift. The LP approach takes advantage of the computer’s computational capacity to solve a set of inequality constraints while minimizing an objective function. The objective function incorporates fuel cost, start-up cost, and wear and tear penalty. The constraints limit the search space to one that satisfies customer demand. In the LP approach, predicting customer demand was not considered. Instead, fuzzy set theory was used to handle uncertainty of customer demand. The system was implemented using the HYPER LINDO package (trademark of Lindo Inc.) and it took significant time to generate a solution. An evaluation version of the software was adopted to carry out the experiments. Due to the limitation on number of constraints supported in the evaluation version, only 6 h per shift were considered in the LP experiments in Sun et al. (2000).
ARTICLE IN PRESS 116
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
2.3. Genetic algorithms GAs are search algorithms based on the mechanics of natural selection and natural genetics. They are guided randomized procedures that efficiently exploit historical information to speculate on the new search points with expected improved performance (Goldberg, 1989). GAs are different from other optimization and search procedures in four ways (Goldberg, 1989): (1) GAs work with a coding of the parameter set, not the parameters themselves. (2) GAs search from a population of points, not a single point, hence different search regions are explored simultaneously. (3) GAs use an objective function to judge quality of solutions instead of using derivatives or other auxiliary knowledge. (4) GAs use probabilistic transition rules instead of deterministic rules. The probabilistic transition rules construct new individuals from the current population based on likelihood that a new population improves upon the performance of the old population. A simple GA is composed of three operators: (1) selection, (2) crossover, and (3) mutation. Selection is a process in which individuals are copied to the mating pool according to their fitness. In this study, this operation was implemented with a Roulette Wheel where individuals with higher fitness occupy more slots. After selection, crossover proceeds to create offsprings. Members of the mating pool are mated at random. This study employed single point crossover and single bit mutation. For each pair of strings, an integer position k is selected uniformly at random and two new strings are created by swapping all characters after position k. Next, mutation is conducted to protect against irrecoverable loss of good genes from the crossover process. It is implemented as an occasional random alternation of the value of a string position. The above-mentioned operator algorithms were chosen because they are the simplest and most commonly used in GA applications. Alternative algorithms will be examined in the future. A discussion on different possible choices of genetic operators is presented below. In tournament selection, fitness values are compared and the best individuals are selected. Unlike the Roulette Wheel approach, this scheme does not reserve potential good genes that are hidden in an undesirable individual. Another approach called Stochastic Remainder uses probabilities to select individuals to be included in the mating pool. Since strong genes usually dominate the population, this approach could give similar result as in the Roulette Wheel approach. However, Stochastic Remainder may have a higher probability of population diversity than the Roulette Wheel approach (Goldberg,
1989). The shape of the problem space often determines which selection algorithm would give the best performance. Two alternative algorithms for the mutation operator are Gaussian and multi-point mutation. Gaussian mutation, which adds a random value from a Gaussian distribution to each element of an individual to create a new offspring, is not applicable for this problem because the individual consists of binary bits. In our opinion, multi-point mutation does not show any obvious advantages over one-point mutation. As the problem under study involves a shift of eight periods each, one way to speed up the search process could be the implementation of an eight-point crossover operator. Since different periods are only loosely related, each period should evolve separately by crossing over with the same period of other individuals. 2.4. Genetic algorithms and constrained optimization In order to solve problems with constraints, an optimization method must handle both feasible and infeasible solutions in a population. Comprehensive surveys of the methods can be found in Dasgupta and Michalewicz (1997) and Coello (2002). Although there are a large number of methods, it is often unclear which method is applicable in a particular situation. In practice, simple penalty functions often prevail in the GA community. Similarly in this study, a static penalty function is used to handle constraints by extending the objective function with penalty terms to penalize infeasible individuals. Consider the following constrained problem as an example: min:
gðxÞ
s: t: where
hi ðxÞX0 i ¼ 1; 2; . . . ; n , x is the controllable vector:
The original problem can be transformed to an unconstrained form as follows: min
gðxÞ þ
n P
ci min2 ð0; hi ðxÞÞ;
i¼1
where
ci is the penalty coefficient:
For each constraint, a penalty term is added to the objective function. If the constraint is satisfied, or hi ðxÞX0, then the term will be evaluated to zero and no penalty is applied. Otherwise, the square of the violating amount is added. The square is used to separate the penalty terms from one another. 2.5. Relevant literature GAs have been applied to a variety of optimization and control problems. Some of their applications in energy engineering are briefly described as follows.
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Goldberg (1989, Chapter 4) had applied GAs to natural gas pipeline control since 1983 but his objectives were different from ours. The objective was to minimize the power consumed subject to maximum and minimum pressure and pressure ratio constraints. Square pressure difference across each compressor station was chosen as the control variable. The constraints were incorporated into the system using a quadratic, external penalty method. The study showed that near optimal results were found quickly after searching an infinitesimal portion of the space. Dasgupta (1997) presented an optimal real-time scheduling system for power generation. The objective of scheduling was to generate power to meet the load with minimum cost of generation. The cost was the sum of start-up, banking, and expected running costs subject to the demand and spinning reverse. The constraints included minimum down time and up time of generating units. A penalty function was added to the objective function to penalize a solution when a constraint is violated. The study suggested that the GA based method is a feasible alternative approach for this type of optimization problem. However, a drawback of this approach is the computational time needed to evaluate the population in each generation. Lu et al. (2005) presented a GA approach to optimize the in-building section of the centralized heating, ventilation and air-conditioning (HVAC) systems that consist of indoor and chilled water loops. The approach aimed to optimize energy consumption of the overall system rather than those of individual components, and unknown parameters were calculated using mathematical models and a neuro-fuzzy inference system. The objective was to minimize the sum of power consumption for chillers, pumps and fans, which was subject to five constraints. Two of these constraints involve structures that make it difficult to obtain their derivatives, rendering it difficult to solve them with traditional techniques. The GAs handled these constraints with penalty functions; and the results showed that the method significantly improved system performance. Osman et al. (2004) proposed a methodology to optimize power flow as a step towards optimal operation and control of power generation. The methodology consists of two GAs and handles constraints by introducing repair functions for infeasible solutions. The first GA searches for an initial reference point that violates fewer constraints than a desired precision. The second one uses the feasible reference point to repair the infeasible solution while searching for an optimal solution. The methodology was applied successfully on a six-bus sample power system. Baran et al. (2005) analyzed six different multiobjective evolutionary algorithms in solving an optimal pump-scheduling problem. The six algorithms are: (1) Strength Pareto Evolutionary Algorithm (SPEA), (2)
117
Non-Dominated Sorting GA (NSGA), (3) Non-Dominated Sorting GA 2 (NSGA2), (4) Controlled Elitist Non-Dominated Sorting GA (CNSGA), (5) Niched Pareto GA (NPGA), and Multiple Objective GA (MOGA). The authors recommended SPEA for the problem. The study differed from Goldberg (1989), Dasgupta (1997), Lu et al. (2005) and Osman et al. (2004) in that it handled multiple objectives separately instead of bundling them into a single objective function. Instead of generating a single solution, the optimization process produced a set of solutions that traded off different objectives including minimizing electric cost, maintenance cost, maximum power peak, and level variation in a reservoir. A repairing algorithm was used to ensure that the solutions generated fall within a feasible range.
3. Neural networks for BHP requirement prediction In order to determine an optimal compressor schedule, one of the factors to consider is customer demand which is difficult to predict. Instead, a heuristic method for calculating the required BHP, which can be used as a proxy for customer demand, was adopted in this study. According to Sun et al., 2000 and Uraikul et al., 2000, the BHP requirement for satisfying customer demand can be estimated using flow rates from the two gas stations with the following heuristic equation: BHP requirement ¼ max ð0; 2:77411 ðSt:Louis Flow þ Melfort FlowÞ 1132Þ.
ð1Þ
Since the independent variables of St. Louis and Melfort flow rates in Eq. (1) are unknown at the beginning of a shift, the neural network technique was adopted to predict flow rates of the current shift based on those of the last shift. Fig. 4 shows the overall structure of the system, which includes two main parts: the neural networks and the GAs. The neural networks take past flow rates from Melfort and St. Louis stations and predict the future flow rates. The outputs of the neural networks are then used for calculating the BHP requirements, which represent the customer demand and in turn serve as part of the input for the GAs. The GAs optimize compressor schedule based on several parameters including BHP requirements, compressor technical parameters, penalty factors and GA search parameters. The compressor technical parameters include fuel cost, start-up cost, and horsepower rating. Penalty factors decide the importance of each type of penalty, including wear and tear penalty, oversupply penalty, and noncovered penalty. These penalties will be discussed further in Section 4.1.2. The GA parameters govern the search process.
ARTICLE IN PRESS 118
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 4. Overall structure of the compressor scheduling system.
Neural network models were trained on historical time-series data of flow-rate recorded at the two stations of St. Louis and Melfort. The data were collected in two periods. The first period was from the cold season of December 3rd, 2001 to March 13th, 2002, and the second period was from the hot season of May 28th to August 8th, 2002. Three alternative neural network structures were considered. In all three alternatives, the Melfort and St. Louis models were developed separately. There are eight input variables in each model, each of which represents the flow rate during the period of an hour in the previous shift. However, the alternative neural network structures differ in the number of output variables and number of neural networks involved. The first structure consists of eight output variables, each of which is a flow rate in the current shift. While this model has the advantage that all the flow rates can be calculated in one step, a weakness is that the connection between the input layer and the hidden layer must account for all eight outputs and hence its underlying structure is more complex. The second structure aims at reducing this complexity by using eight separate neural networks; each neural network predicts the flow rate for one of the 8 h in a shift. The disadvantage of this approach is the large number of neural networks involved; altogether, 32 analytical and 16 final models must be built. In the third alternative structure, only the model that predicts 1 h ahead is necessary and this model will be used recursively to predict 8 h ahead. This type of model has low complexity and can adequately explain training data. However, since the models are trained for short term and used for long term forecast, the generalization may be unreliable. In selecting among the three alternative structures, the second alternative was rejected because it involved implementing and managing a large number of neural
Table 1 Models’ validation root mean square errors Station
Training season
Validation season
Number of outputs
RMSE
Melfort Melfort Melfort Melfort St. Louis St. Louis St. Louis St. Louis
Cold Cold Hot Hot Cold Cold Hot Hot
Hot Hot Cold Cold Hot Hot Cold Cold
1 8 1 8 1 8 1 8
113.55 55.93 91.43 74.05 101.53 118.35 225.10 83.92
networks. The remaining two alternatives were implemented and analyzed. Since characteristics of the hot vs. cold seasons are significantly different, the data sets for one season were used to develop models while the other sets were used to validate the models’ generalization. As can be seen in Table 1, the first alternative (highlighted rows) generated lower root mean square errors (RMSE) on the validation data sets for 8 h ahead, and hence it was used for developing the final models. The data sets for both seasons were then combined to build the final models which were incorporated into the GA system. The training results are summarized in Table 2. During training, the data was re-sampled in different ways into two portions for training and testing; and the testing portion ensured not only that the model fitted well to the data but also that it could be generalized to cover unseen data. Both training and testing errors and correlative coefficients were recorded, and plots of model-generated and observed data points for 1 hand 8 h ahead are displayed in Figs. 5 and 6, respectively. As expected, prediction accuracy decreased as the terms of prediction increased. For example, in Table 2, the testing RMSE of the St. Louis model
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
119
Table 2 Training statistic of final models Hours ahead
Melfort
St. Louis
Training
8 7 6 5 4 3 2 1
Testing
Training
Testing
RMSE
CORR_COEFF
RMSE
CORR_COEFF
RMSE
CORR_COEF
RMSE
CORR_COEF
50.88 50.84 49.01 46.87 54.16 42.61 38.43 28.94
0.8392 0.8442 0.8545 0.8648 0.8778 0.8909 0.9115 0.9501
64.79 61.51 59.63 57.14 53.07 48.75 43.31 34.06
0.7485 0.7742 0.7875 0.8033 0.8333 0.8603 0.8906 0.9323
62.97 58.81 57.87 56.43 52.26 47.75 41.52 31.43
0.9152 0.9292 0.9315 0.9343 0.943 0.9526 0.9642 0.9796
70.02 64.25 62.49 60.39 57.22 52.35 45.87 34.06
0.9211 0.9241 0.9284 0.933 0.9398 0.9497 0.9615 0.9785
increased from 34.06 to 70.02 m3/h when the term of prediction increased from 1 to 8 h ahead. This can also be observed from Figs. 5 and 6. The dots in Fig. 5 gather closely to the line while the dots in Fig. 6 are more scattered. The training and testing correlative coefficients (CORR_COEFF in Table 2) of the St. Louis model was higher than those of the Melfort models, which means that the St. Louis model both explains the training data better and has higher generalization ability than the Melfort model. On the other hand, the lowest testing correlative coefficient of the Melfort model was 0.7485, which indicates that the trained model can explain at least 74.85% of the unseen data. This is an acceptable result. The neural network models were implemented in Neur-Online (NOL) Studio (trademark of Gensym Corp., USA) and were imported to the main program as ActiveX components. From the neural network modeling, the flow rates and hence the BHP requirements are obtained, which becomes input to the second subproblem of compressor scheduling. Next, we present representation of the problem of compressor scheduling for GA modeling.
4. Genetic algorithms for optimization of gas compressor schedule 4.1. Problem representation 4.1.1. Chromosome The first step in developing a GA is to encode the problem space of operational schedules. The entire operational schedule was encoded in a single chromosome. There are five controllable compressors operating in a time span of 8 h. The change of compressor states is assumed to occur only at the beginning of each hour. Hence, there are a total of 5 times 8 or 40 bits in a chromosome string. Each bit will have a value of 1 if the compressor is turned on and 0 if the compressor is
turned off. Since there are two choices of value for the first bit, and for each choice of the first bit there are two choices of the second bit, and so on, the total number of possible solutions to this optimization problem is 24041012. It could be difficult for a traditional method to explore such a large search space; hence, GAs, being an efficient heuristic search technique, are employed to analyze the problem. Two alternatives organizations were considered for the order of the chromosome. A chromosome can either contain eight sub-strings for eight periods of five compressors’ states, or five sub-strings of compressors’ states in eight periods. In the first alternative, if a chromosome is represented with an array of bits x[1..40], x[1] to x[5] will represent the states of the five compressors in the first hour, x[6] to x[10] will represent the states of the same compressors in the second hour, and so on. In each sub-string, the first three compressors run on natural gas and the last two run on electricity. In the second alternative, x[1] to x[8] will represent the states of the first compressor in the eights periods of an hour each, x[9] to x[16] will represent the states of the second compressor in the eight periods, and so on. The order difference does not affect fitness of a chromosome; however, it could make a difference in the search since it is more difficult to break up genes that are close together in a crossover process. Since related genes should be placed close to one another, the relationships among genes are examined. Referring to the cost and penalty functions in the next section, a sequence of states of a compressor in a shift influences start-up cost and wear-and-tear penalty. However, calculation of start-up cost relies only on two consecutive periods, and wear-and-tear penalty involves five periods (the calculation of wear-and-tear penalty will be explained in Section 4.1.2). A combination of compressors in a period, on the other hand, influences penalties for oversupply and non-covered demand in that period. Since the non-covered demand penalty is a hard constraint that should not be violated,
ARTICLE IN PRESS 120
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 5. St. Louis Model—predicted vs. observed—1 h ahead.
it is more important than the other costs and penalties. Hence, it is concluded that the relationship between compressors in a period is more important than the sequence of a compressor’s states in eight periods, and the first alternative organization of chromosome order is preferred. The binary coding reduces the number of constraints from seven in the LP representation approach (Sun et al., 2000) to one in the GA approach. The reason for this reduction is explained as follows. In the pipeline problem, compressor states have discrete values of 0 or 1 that can be easily represented as binary bits in a chromosome. However, since discrete parameter values cannot be directly represented in the LP approach, a number of extra inequality constraints were added. Without these extra constraints, the LP system would generate infeasible solutions. For example, in Sun et al.
(2000), to represent the binary state of natural gas compressor i in period p (denoted as Ni, p), a linearization variable xi;p ¼ N i;p N i;p1 was introduced and the following three extra constraints were added. N i;p þ N i;p1 xi;p p1 8i;p ; N i;p N i;p1 þ 2xi;p p0 8i;p ; N i;p þ N i;p1 xi;p X0 8i;p : These three extra constraints are always true in the GA setting and therefore unnecessary. 4.1.2. Fitness evaluation The objective function is specified as the sum of fuel cost, start-up cost, wear and tear penalty (i.e. cost for not letting a compressor run for at least three consecutive hours), oversupply penalty (i.e. cost for
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
121
Fig. 6. St. Louis Model—predicted vs. observed—8 h ahead.
oversupply of gas), and non-covered penalty (i.e. cost for not satisfying customer demand). Since individuals with higher fitness are preferred, the fitness function is expressed as the negation of the objective function.
s(i) d(p)
fitnessðxÞ ¼ ½Fuel_CostðxÞ þ Startup_CostðxÞ
c1(i)
þ Wear_Tear_PenaltyðxÞ þ Oversupply_PenaltyðxÞ þ Non_Covered_PenaltyðxÞ. Let us denote the problem parameters as follows: x[1..40] Array representing a chromosome, where each element x[j], jA[1, 40] represents the state of the (j mod 5)th compressor in the (floor(j/5))th period f(i) One-hour fuel cost for compressor i, iA[0, 4]
h(i)
c2 c3
Start-up cost for compressor i, iA[0, 4] BHP requirement in period p (Proxy for customer demand), pA[1, 8] Horsepower rating of the ith compressor, iA[0, 4] Wear and tear penalty factor for compressor i, iA[0, 4] Oversupply penalty factor Non-covered penalty factor
The partial costs and penalties are calculated as follows:
Fuel cost: The fuel cost is calculated as total fuel cost consumed by all compressors of the two gas stations during a shift. The fuel consumption rates of
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
122
compressors are fixed and different. The cost for running a natural gas compressor is higher than that for an electrical compressor. Fuel_CostðxÞ ¼
40 X
The pseudo-code can be translated into the following equivalent equation: Wear_Tear_PenaltyðxÞ ¼
f ðj mod 5Þx½j.
Start-up cost: The start-up cost is accumulated each time a compressor is turned on from its idle state. Start-up cost for an electrical compressor is higher than that for the natural gas compressors. Start_CostðxÞ ¼
40 X
sðj mod 5Þ maxðx½j x½j 5; 0Þ.
i¼6
Wear and tear: The combination of compressors should be selected in such a way that the number of switching is minimized. After a compressor is turned on, it should be kept in operation for as long as possible. Similarly, when a compressor is shut off, it should be kept inactive for as long as possible. In this study, the minimum operative and inoperative period is set to be 3 h.
c1 ½j mod 5
j¼6
x½j x½j þ 5 þ x½j x½j þ 5 ! 0 mþ2 Y X 2 ð x½j þ 5 kÞ mod 3 .
j¼1
35 X
m¼2 k¼m
This equation is non-linear and cannot be solved using LP. Sun et al. (2000) provided a linear equation to describe the wear and tear factor, which punishes any three consecutive sequences that have a sum of less than three. However, the equation wrongly punishes the free flow solution in which no compressors are running. Another flaw in the approach reported in Sun et al. (2000) is that the wear and tear penalty did not take past states of the compressors into consideration.
Oversupply penalty: This penalty applies when there is surplus gas delivered to the customers. The oversupply amount of gas will be wasted and is assigned a penalty. OverSupply_PenaltyðxÞ
The wear and tear penalty is implemented as follows. The penalty for each compressor is calculated by examining the center of a moving window of five successive periods w[1..5]. If the center w[3] is the same as either of its two predecessors w[2] and w[1], successors w[4] and w[5], or immediate neighbors w[2] and w[4], it should not be punished. Otherwise, the penalty will be calculated based on the center’s immediate neighbors. If we do not have enough information to determine whether or not we should penalize a period, we will not punish the solution. For example, if information on the previous shift is not available then the first period of the current shift should never get penalized for wear and tear. Let us denote st[1..8] as the array of the states of a compressor during a shift. In the current implementation, the eight periods 1–8 are extended to 10 periods 0–9 where states st[0] ¼ st[1] and st[9] ¼ st[8]. The pseudo-code for calculating the wear and tear penalty is as follows, where p refers to index in [1..8].
For period p from 1 to 7 If (not ((st[p] ¼ st[p1] and st[p] ¼ st[p2]) or (st[p] ¼ st[p1] and st[p] ¼ st[p+1]) or (st[p] ¼ st[p+1] and st[p] ¼ st[p+2]))) then penalty ¼ penalty+c1*(|st[p]st[p1]| +|st[p]st[p+1]|) End if End for
¼ c2 max 0;
40 X
!!2 hðj mod 5Þx½j dðpÞ
.
j¼1
Penalty for non-covered demand: This penalty applies when the gas supply cannot satisfy customer demand. If the supply is sufficient, no penalty is applicable; otherwise, the part of demand not covered by supply will be punished. The penalty factor c2 of this penalty should be much higher than the penalty factor c1 of the oversupply penalty because if demand is not covered, customer satisfaction cannot be ensured. Non_Covered_PenaltyðxÞ ¼ c2 max 0; dðpÞ
40 X
!!2 hðj mod 5Þx½j
.
j¼1
The constraints in the GA representation are more compact than those in the LP approach. In the GA approach, all constraints are bundled in one value while there are eight inequalities for the eight periods in the LP approach. 4.1.3. Genetic algorithm process The simple GA that was introduced briefly in Section 2.3 was used to solve the compressor-scheduling problem. The algorithm started with a small number of random individuals in the initial population. At each
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
generation, the better individuals were chosen with higher probabilities for creating successors. The best individuals from each generation were kept. The runtime of the GA for 1000 generations took only a few minutes. This efficiency in runtime can possibly be explained as follows: (1) the GA approach only examines a small number of individuals each time, (2) it only searches in some parts of the search space that are deemed to be closer to the global optima instead of searching the entire problem space as in exhaustive search, and (3) the time for calculation and comparison of fitness function is likely to be less than the time needed to solve equations.
a value of 30 or 50, the number of generations varied from 100 to 1000, the crossover rate had values of 0.6, 0.8 or 1.0 and the mutation rate had values of 0.0005, 0.004, and 0.03. The GAs were executed from three to five times with each combination of these values. Although the results are different in each run due to randomness in the initial population, the following observations can be made.
4.2. Discussions This section discusses the GA search process from two pespectives. First, the impact of parameter values on the search process of the GA will be examined. Secondly, the penalty factor will be discussed. The GAs were executed with the fixed parameters of fuel cost, start-up cost, and horsepower rating. The values for the parameters are listed in Table 3. Other parameters may vary. As can be seen from Table 3, the two electrical compressors (E1 and E2) have the same technical specification. Similarly, so do the two natural gas compressors (N2 and N3) at Melfort station. The natural gas compressors have higher horsepower rating but also require more fuel to operate. 4.2.1. Genetic algorithm parameters In order to examine the impact of parameter values on the search process, the parameter values in the GA were varied as follows. The population size was assigned
Table 3 Fixed parameters Compressor
Fuel cost
Start-up cost
Horsepower rating
N1 N2 N3 E1 E2
400 350 350 100 100
120 120 120 200 200
650 600 600 250 250
123
The change in population size and crossover rate did not affect the genetic search in any way. When the BHP requirements were small, the search converged quite quickly. However, when the BHP requirements increased, the best solution usually appeared later. Hence, a large maximum number of generations are always recommended. The mutation process is needed to recover useful genetic material that was accidentally lost over generations. Hence, as can be observed in Fig. 7, when the mutation rate was too small (0.0005), the search was trapped and could not return to its initial good solutions in later generations. The initial best solution has the fitness values in the range of 80–90 but since the fitness drop to 50 and under in a bad transition, it cannot regain its initial values. On the other hand, if the mutation rate is too high (0.03) as in Fig. 8, the search oscillated randomly. Finally, when the rate is right, the fitness improved through generations (see Fig. 9). The charts on the left of Figs. 7–9 represent the scaled fitness of the best solutions through generations and the charts on the right represent the average solutions.
4.2.2. Penalty factors The penalty factors or coefficients have significant impact on the quality of solutions. Among the factors, the penalty factor for non-covered demand is the most important, because it ensures customer satisfaction and represents a hard constraint that should not be violated in the final best solution. While some infeasible solutions in the population can provide diversity, the number of infeasible solutions should not be too large.
Fig. 7. Effect of a small mutation rate.
ARTICLE IN PRESS 124
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 8. Effect of a large mutation rate.
Fig. 9. Effect of a good mutation rate.
It is observed that the penalty factor for non-covered demand (c3) should be proportional to the BHP requirements. If the BHP requirements are small, the feasible search space is large, so even if the penalty is small, there are still some good feasible solutions in the population. However, when the BHP requirements are large, the penalty should be stricter in order to maintain feasibility of population. Our experiment showed that when the BHP requirements are 2000 or above, the penalty factor for noncovered demand must be 3,000,000 or above in order to have a feasible solution. However, even with such a high factor, the number of feasible solutions in all generations was low. Nevertheless, the average number of constraint violations of each chromosome seemed to reduce from 4 to 6 each at early generations to 1–3 in later generations. The disadvantage of a high non-covered penalty factor is that it will overshadow the importance of the wear and tear and other cost and penalties. In the case where there are mixed low and high BHP requirements (e.g. low in periods 1–4 and high in periods 5–8), this problem can produce costly solutions. For example, let us look at a solution when the BHP requirement is high as shown in Table 4. The highest BHP requirement in the eight periods is 2000, which means the non-covered penalty factor is high. This overshadows for example, the wear and tear cost. Since E1 and E2 are of the same capacity, the solution is not optimal because of the unnecessary back and forth switches between E1 and E2 from period 3 to 4 and from period 5 to 6. The solution
Table 4 A genetic algorithm solution Period
BHP requirement
N1
1 2 3 4 5 6 7 8
800 1000 1100 1500 1600 2000 1300 1100
On
On On On On On
N2
On On On On On
N3
E1
E2
On On On
On
On On On
On On On On
On On On On
On On
is better if E1 remains off and E2 remains on in periods 4 and 5. Hence, each period should have its own penalty factor for non-covered demand instead of sharing the same penalty factor as in the current implemented system. This can be considered in future work.
5. Results: comparison with past solutions In this section, the solutions from the three methods of expert system, LP, and GAs will be comprared in two different cases with different sets of BHP requirements. As all the solutions satisfy the hard constraint of ensuring customer satisfaction, the comparison was made in terms of wasted energy and compressor switches between consecutive periods. As mentioned in Section 2.1, the evaluation version of the HYPER LINDO software package for LP could
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
handle a maximum of only 200 constraints, and only six periods can be optimized at each run. To make a reasonable comparison, the BHP requirements from Sun et al. (2000) were used for periods 1–6 and two random BHP requirements were used for periods 7 and 8. Table 5 displays the solutions from the three methods of expert system, LP, and GAs. The GA solution wasted less energy than the expert system but is not as efficient as the LP solution. For example, in the fifth period, the expert system used all five compressors while the GA turned off an electrical compressor. The LP approach went a step further when it substituted a gas compressor in the GA solution with a less powerful electrical compressor. In the first six periods, the total horsepower supplied with the schedules generated by expert system, LP and GA are 10950, 8700, and 9750, respectively, while the necessary horsepower was only 8000. Given the discussion in Section 4.2.2, the wasted resources can be accounted for as follows: since the BHP requirement is high, the non-covered demand penalty factor was too high and hence, other costs and penalties were neglected. However, in terms of wear and tear penalty, the expert system provided the most efficient solution because in the first six periods, the BHP requirements increase monotonously and hence compressors were only added and not turned off. The only penalty that the expert system was subjected to is when compressor N3
125
was turned off in period 7 after only two running periods. The GA solution was inefficient in terms of wear and tear cost, because compressor switches were frequent in almost every period of operation. On the other hand, the LP seemed to ignore the wear and tear factor as they allowed the OR operator in its set of solutions. As an example, if the combinations of the first two periods are {E1, N2} and {E1, E2, N2}, then there are no switches; but if they are {E1, N2} and {E1, E2, N3}, then there is a switch between two equivalent natural gas compressors N2 and N3. Another experiment was conducted; this time the BHP requirements were predicted based on previous flow using neural networks. No results from LP are available for this set of BHP requirements. The results from GA are shown in Table 6. Again in this case, we observe that the solution from GA wasted less resource than that from the expert system. The expert system experienced a sudden shift in the amount of energy generated when the BHP requirement crossed the 800 level. In periods 1 and 2, when the BHP requirements were close but still less than the upper bound of 800 HP, the expert system was very efficient. However, it was least effective when it barely crossed the boundary of 800 HP because the total energy generated by the expert system approach was 1500 HP, which means (1500–800) or 700 units of horsepower was
Table 5 Comparison of the three approaches Period
BHP requirement
Expert system
Linear programming
Genetic algorithm
1 2 3 4 5 6 7 8
800 1000 1100 1500 1600 2000 1300 1100
E2, E2, E2, E2, E2, E2, E2, E2,
(E1|E2), (N2|N3) E1, E2, (N2|N3) E1, E2, (N1|N2|N3) E1, E2, N1, (N2|N3) E1, E2, N1, (N2|N3) (E1|E2), N1, N2, N3
N1, N3, N2, N1, N1, N1, N1, N1,
N1, N1, N1, N1, N1, N1, N1, N1,
N2 N2 N2 N2, N2, N2, N2, N2
E1 N3, E1 N3, E1 E1
N3, E2 E1, E2 N3, E2 N2, E1 N2, N3, E1 N2, N3, E2 N3, E1, E2 N2, N3, E1
Table 6 Genetic algorithm solution in the case of medium BHP requirement Period BHP requirements Expert system solution Genetic algorithm solution Total energy generated by ES Total energy generated by GA 1 2 3 4 5 6 7 8 Total
796.09 799.73 803.98 801.71 811.94 813.71 818.88 826.88 6472.92
E2, E2, E2, E2, E2, E2, E2, E2,
N1 N1 N1, N1, N1, N1, N1, N1,
N2 N2 N2 N2 N2 N2
N1, N1, N2, N1, N3, N2, N2, N3,
E1, E2 N3 E1 N3, E1 E2 N3 E1 E1, E2
900 900 1500 1500 1500 1500 1500 1500 10800
1150 1250 850 1500 900 1200 850 1100 8800
ARTICLE IN PRESS 126
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
wasted. In terms of wear and tear cost, the GA solution was also inefficient in this experiment because of frequent switching.
population replacement or De Jong’s crowding model (Goldberg, 1989).
Acknowledgment 6. Conclusions and future work The genetic algorithms (GAs) in this study provided an efficient way to search for a compressor operation schedule. The running time for 1000 generations of 30 individuals each is less than 5 min which is significantly shorter than that of the linear programming approach. However, it is difficult to determine if a near optimal solution has been found or whether the search has been trapped in a local minimum. The quality of the solution also depends heavily on the initial random population. The fact that GAs use payoff information, not derivatives or other auxiliary knowledge, gives it more flexibility than other traditional methods. GAs can handle non-linearity in the same way as linearity. However, this heuristic also involves a weakness in that some information of the function to be optimized is lost in the process of calculating payoffs. The approach to handle constraints in this study is not yet efficient. It depends largely on a number of penalty coefficients that need to be tuned. In future research, the use of repairing technique to repair infeasible solutions will be considered. A simple repairing algorithm could be as follows. Until the demand is met, keep turning on more compressors, which are selected at random. In the case of small feasible search space such as when high BHP level is required, the repairing approach is especially helpful since it ensures that a feasible solution is found. Another possible direction for future work is to add other constraints such as maximum operating hours or real-time failure. The genetic algorithms can also be extended with more complex features such as elitist
The authors would like to acknowledge the generous support of a Strategic Grant from the Natural Sciences and Engineering Research Council of Canada. References Baran, B., von Lucken, C., Sotelo, A., 2005. Multi-objective pump scheduling optimization using evolutionary strategies. Advances in Engineering Software 36, 39–47. Coello, C.A., 2002. Theoretical and numerical constraint-handling techniques used with evolutionary algorithms: a survey of the state of the art. Computer Methods in Applied Mechanics and Engineering 191 (11–12), 1245–1287. Dasgupta, D., 1997. Optimal Scheduling of thermal power generation using evolutionary algorithms. In: Dasgupta, D., Michalewicz, Z. (Eds.), Evolutionary Algorithms in Engineering Applications. Springer, Berlin, Heidelberg, pp. 317–328. Dasgupta, D., Michalewicz, Z., 1997. Evolutionary algorithms—an overview. In: Dasgupta, D., Michalewicz, Z. (Eds.), Evolutionary Algorithms in Engineering Applications. Springer, Berlin, Heidelberg, pp. 3–28. Goldberg, D.E., 1989. Genetic Algorithms in Search, Optimization, and Machine Learning. Addition-Wesley Publishing Company Inc., Reading, MA. Lu, L., Cai, W., Xie, L., Li, S., Soh, Y.C., 2005. HVAC System optimization—in-building section. Energy and Buildings 37, 11–22. Osman, M.S., Abo-Sinna, M.A., Mousa, A.A., 2004. A solution to the optimal power flow using genetic algorithm. Applied Mathematics and Computation 155, 391–405. Sun, C.K., Uraikul, V., Chan, C.W., Tontiwachwuthikul, P., 2000. An integrated expert system/operation research approach for the optimization of natural gas pipeline operations. Engineering Applications of Artificial Intelligence 13, 465–475. Uraikul, V., Chan, C.W., Tontiwachwuthikul, P., 2000. Development of an expert system for optimizing natural gas pipeline operations. Expert Systems with Applications 18, 271–282.
Engineering Applications of Artificial Intelligence 19 (2006) 113–126 www.elsevier.com/locate/engappai
Applications of artificial intelligence for optimization of compressor scheduling Hanh H. Nguyen, Christine W. Chan Faculty of Engineering, University of Regina, Regina, Saskatchewan, Canada S4S 0A2 Received 2 April 2005; received in revised form 13 June 2005; accepted 16 June 2005 Available online 10 October 2005
Abstract This paper presents a feasibility study of evolutionary scheduling for gas pipeline operations. The problem is complex because of several constraints that must be taken into consideration during the optimization process. The objective of gas pipeline operations is to transfer sufficient gas from gas stations to consumers so as to satisfy customer demand with minimum costs. The scheduling involves selection of a set of compressors to operate during a shift. The scheduling decision has to be made so as to satisfy the dual objectives of minimizing the sum of fuel cost, start-up cost, the cost of gas wasted due to oversupply, and satisfying minimal operative and inoperative time of the compressors. The problem was decomposed into the two subproblems of gas load forecast and selection of compressors. Neural networks were used for forecasting the load; and genetic algorithms were used to search for a near optimal combination of compressors. The study was conducted on a subsystem of the pipeline network located in southeastern Saskatchewan, Canada. The results are compared with the solutions generated by an expert system and a fuzzy linear programming model. r 2005 Elsevier Ltd. All rights reserved. Keywords: Genetic algorithms; Neural networks; Compressor scheduling
1. Introduction Evolutionary scheduling for gas pipeline compressor is a complex problem because of constraints that must be taken into consideration during the optimization process. Such constraints include operation cost, maintenance cost, and customer satisfaction. This paper reports on a feasibility study on applying Artificial Intelligence technologies to natural gas pipeline operations. The focus of the study was on the St. Louis East gas pipeline distribution system, a subsystem of the pipeline network located in southeastern Saskatchewan. The system consists of two stations located at Melfort and St. Louis. The two stations use compressors to Corresponding author. Tel.: +1 306 585 5225; fax: +1 306 585 4855. E-mail addresses: [email protected] (H.H. Nguyen), [email protected] (C.W. Chan).
0952-1976/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.engappai.2005.06.008
generate sufficient pressure to transport gas to the surrounding consumption areas of St. Brieux, Nipawin and Hudson Bay. The dispatchers’ task is to adjust the compression level in order to generate the necessary pressure while not wasting energy. When the customer demand increases, a dispatcher adds compression to the pipeline system by turning on one or more compressors. On the other hand, the dispatcher turns off one or more compressors to reduce compression in the pipeline system when the customer demand decreases. This decision has a significant impact on effectiveness of the natural gas pipeline operation as well as on operating costs. The objective of this study is to use automation techniques to aid the dispatcher in optimizing natural gas pipeline operations in order to satisfy customer demand with minimal operating costs. A dispatcher needs to know in advance when a large volume of gas is required and what way is the optimal to ensure it is
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
available. Previous studies were conducted in which expert system and linear programming (LP) technologies were employed to solve the problem. Here, we present a genetic algorithm (GA) system that generates workable schedules of compressor operation based on the current situation. Since the problem can be decomposed into two sub-problems of demand forecast and compressor selection, neural networks are employed to tackle the former problem of demand forecast and GAs are used for addressing the problem of compressor optimization. This paper is organized as follows. Section 2 provides some background on the past solutions for the problem and the GA technique in general. Section 3 describes the first subtask, as well as development and results of demand prediction using neural networks. Section 4 discusses the second subtask, which involves analysis and implementation of GAs for optimization of compressor scheduling. Section 5 presents comparison of results generated by the three methods. Section 6 gives some conclusion and discusses some possible future work.
2. Background 2.1. Problem domain This project was part of a study on automation of natural gas pipeline operations conducted jointly with a gas pipeline company in Saskatchewan, Canada. The objective of the study was on optimization of compressor scheduling and the study focused on a section of the Saskatchewan gas pipeline around the St. Louis East compressor station in the province of Saskatchewan, Canada. A schematic of the St. Louis East system is shown in Fig. 1. The system consists of two compressor stations, Melfort and St. Louis. These compressor stations supply natural gas to two customer locations, Nipawin and Hudson Bay. In St. Louis, there are three compressor units; two of them are electrical compressor units and
Melfort compressor station
Nipawin
Hudson Bay Electricity compressor Gas compressor
Fig. 1. Schematic of the St. Louis East system.
Industrial customer Time
2
Pressure
1
Pressure
St. Louis compressor station
the other is a gas compressor unit. In Melfort, there are two gas compressors. An electrical compressor unit provides 250 hp and a gas unit provides 600 or 650 hp. The demand for natural gas from customers fluctuates depending on the season. In the winter, the demand for natural gas is usually higher than in the summer. In addition, the demand for natural gas also changes depending on the time of day. For example, in the morning, the demand is higher because the customer begins to use natural gas to warm the premise. In the afternoon, the demand is low since the facilities are already heated up in the morning. The customers can also be grouped into three types: industrial, dehydrator, and heat sensitive customers. Each type of customer reflects a different pattern of natural gas consumption as illustrated in Fig. 2. The industrial customer has the same rate of natural gas consumption any time of the day. The dehydrator customer demands two set amounts of natural gas at different periods of time. For example, between 8:00 and 10:00 a.m., the demand is 240 103 m3/day while between 10:00 am to 12:00 am, the demand is 200 103 m3/day. The demand of the heat sensitive customer fluctuates over time and is difficult to predict. In practice, the demand for natural gas can fluctuate from 200 103 m3/day to over 560 103 m3/day within 1 or 2 h. A demand of 200 103 m3/day can be handled from St. Louis with one compressor unit. When the demand exceeds 560 103 m3/day, all units at St. Louis and both units at Melfort are needed. To ensure customer satisfaction, the operator needs to know beforehand and be ready for the largest volume requirement of gas. Otherwise, the system pressures at Nipawin and Hudson Bay will be below the required minimum. The graphs showing supply and demand of gas in the pipelines are generated during operations by the simulator program of the gas pipeline transmission company. The display includes a compressor discharge pressure curve vs. customer station pressure curve as shown in Fig. 3. The two curves in Fig. 3 provide pressure information to dispatchers, and indicate the relationship between demand and supply. The gap between the two curves is called the comfort zone (CZ). If this CZ is wide, then customer satisfaction is
Dehydrator customer Time
3
Pressure
114
Heat sensitive customer Time
Fig. 2. Different natural gas consumption patterns.
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126 6500 6000
Pressure (kpa)
5500 5000 4500 4000 3500 3000
Compressor Discharge Pressure (supply) Customer End-point Delivery Pressure (demand)
2500 0
200
400
600
800
Time (hr)
Fig. 3. Comfort zone.
guaranteed but the cost is high; if it is too narrow, then the operation is cost effective but future customer demand may not be satisfied. Hence, the dispatcher needs to monitor the volume of gas supplied to the customer and maintain a suitable width of the CZ between the two pressure curves. Optimizing compressor scheduling and deciding which of the five compressors to turn on/off is a primary step for increasing pressure in the pipelines while maintaining the optimum operating pressure in the pipeline system. A novice dispatcher may make unnecessary start/stop decisions due to lack of experience. A minimum of 4 months is required to train new dispatchers to operate the system. However, it takes a number of years for them to gain proficiency in operating the system smoothly and cost effectively. 2.2. Compressor scheduling: expert system and linear programming approaches The problem of compressor schedule optimization has been previously tackled with the technologies of expert system (Uraikul et al., 2000) and fuzzy LP (Sun et al., 2000). The expert system tackled the three subtasks of (1) determining the state of the current line pack to decide if compression needs to be added or reduced, (2) determining the required break horsepower (BHP), and (3) selecting the combination of compressors from the two stations to meet the BHP requirement with minimal cost. Since the main emphasis of the expert system was to support the user in the first two tasks of determining the control action and BHP, a rather simplistic approach was adopted for the third task of compressor selection. The system only took into consideration the BHP requirement as the basis for scheduling compressors. The system recommended five combinations of compressors to be turned on in the following order: (1) Free flow: St. Louis E2, (2) 0oBHPo800: St. Louis E2 and N1, (3) 800pBHP o1200: St. Louis E2 and N1,
115
and Melfort N2, (4) 1200pBHP o1600: St. Louis E2, E1 and N1, Melfort N2, and (5) 1600pBHPo2200: St. Louis E2, E1, N1, and Melfort N2 and N3, where St. Louis E1 and E2 are electrical compressors, and St. Louis N1, Melfort N2 and N3 are natural gas compressors. This suggestion is made based on heuristic knowledge of senior pipeline operators. The expert system has been demonstrated using simulation data for various operational scenarios to a panel of expert operators, and they have confirmed that its recommendations were valid. However, no attempts have been made to quantify the performance in terms of component costs. While customer demand satisfaction was an implicit assumption in the experts’ formulation of the rules in the expert system, other component costs were not explicitly considered. As a result, neither the rules nor the validation results reveal how well the expert system handles the other costs. A closer examination on these rules reveals some interesting hidden heuristics. For example, the rules involve only either addition or substraction of compressors, not both at the same time. Given the fact that it is unlikely for BHP requirements to oscillate between two ranges of BHP in a short time period of 3 h or less, the chance of incurring wear and tear on the units is limited. In terms of cost, a disadvantage of the expert system approach is that it enforces safety by sacrificing energy efficiency. The expert system rules are constructed in a way that ensures the upper bound of every range of BHP is always satisfied. Therefore, the rule on compressor selection for that range is not efficient when the actual BHP requirement is in fact closer to the lower bound of the range. A weakness of the expert system is that its decisions are based solely on data from instaneous snap-shots of the pipeline and hence ignore some important information such as demand trend. On the other hand, the LP and GA approaches take into consideration data from the entire duration of a shift. The LP approach takes advantage of the computer’s computational capacity to solve a set of inequality constraints while minimizing an objective function. The objective function incorporates fuel cost, start-up cost, and wear and tear penalty. The constraints limit the search space to one that satisfies customer demand. In the LP approach, predicting customer demand was not considered. Instead, fuzzy set theory was used to handle uncertainty of customer demand. The system was implemented using the HYPER LINDO package (trademark of Lindo Inc.) and it took significant time to generate a solution. An evaluation version of the software was adopted to carry out the experiments. Due to the limitation on number of constraints supported in the evaluation version, only 6 h per shift were considered in the LP experiments in Sun et al. (2000).
ARTICLE IN PRESS 116
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
2.3. Genetic algorithms GAs are search algorithms based on the mechanics of natural selection and natural genetics. They are guided randomized procedures that efficiently exploit historical information to speculate on the new search points with expected improved performance (Goldberg, 1989). GAs are different from other optimization and search procedures in four ways (Goldberg, 1989): (1) GAs work with a coding of the parameter set, not the parameters themselves. (2) GAs search from a population of points, not a single point, hence different search regions are explored simultaneously. (3) GAs use an objective function to judge quality of solutions instead of using derivatives or other auxiliary knowledge. (4) GAs use probabilistic transition rules instead of deterministic rules. The probabilistic transition rules construct new individuals from the current population based on likelihood that a new population improves upon the performance of the old population. A simple GA is composed of three operators: (1) selection, (2) crossover, and (3) mutation. Selection is a process in which individuals are copied to the mating pool according to their fitness. In this study, this operation was implemented with a Roulette Wheel where individuals with higher fitness occupy more slots. After selection, crossover proceeds to create offsprings. Members of the mating pool are mated at random. This study employed single point crossover and single bit mutation. For each pair of strings, an integer position k is selected uniformly at random and two new strings are created by swapping all characters after position k. Next, mutation is conducted to protect against irrecoverable loss of good genes from the crossover process. It is implemented as an occasional random alternation of the value of a string position. The above-mentioned operator algorithms were chosen because they are the simplest and most commonly used in GA applications. Alternative algorithms will be examined in the future. A discussion on different possible choices of genetic operators is presented below. In tournament selection, fitness values are compared and the best individuals are selected. Unlike the Roulette Wheel approach, this scheme does not reserve potential good genes that are hidden in an undesirable individual. Another approach called Stochastic Remainder uses probabilities to select individuals to be included in the mating pool. Since strong genes usually dominate the population, this approach could give similar result as in the Roulette Wheel approach. However, Stochastic Remainder may have a higher probability of population diversity than the Roulette Wheel approach (Goldberg,
1989). The shape of the problem space often determines which selection algorithm would give the best performance. Two alternative algorithms for the mutation operator are Gaussian and multi-point mutation. Gaussian mutation, which adds a random value from a Gaussian distribution to each element of an individual to create a new offspring, is not applicable for this problem because the individual consists of binary bits. In our opinion, multi-point mutation does not show any obvious advantages over one-point mutation. As the problem under study involves a shift of eight periods each, one way to speed up the search process could be the implementation of an eight-point crossover operator. Since different periods are only loosely related, each period should evolve separately by crossing over with the same period of other individuals. 2.4. Genetic algorithms and constrained optimization In order to solve problems with constraints, an optimization method must handle both feasible and infeasible solutions in a population. Comprehensive surveys of the methods can be found in Dasgupta and Michalewicz (1997) and Coello (2002). Although there are a large number of methods, it is often unclear which method is applicable in a particular situation. In practice, simple penalty functions often prevail in the GA community. Similarly in this study, a static penalty function is used to handle constraints by extending the objective function with penalty terms to penalize infeasible individuals. Consider the following constrained problem as an example: min:
gðxÞ
s: t: where
hi ðxÞX0 i ¼ 1; 2; . . . ; n , x is the controllable vector:
The original problem can be transformed to an unconstrained form as follows: min
gðxÞ þ
n P
ci min2 ð0; hi ðxÞÞ;
i¼1
where
ci is the penalty coefficient:
For each constraint, a penalty term is added to the objective function. If the constraint is satisfied, or hi ðxÞX0, then the term will be evaluated to zero and no penalty is applied. Otherwise, the square of the violating amount is added. The square is used to separate the penalty terms from one another. 2.5. Relevant literature GAs have been applied to a variety of optimization and control problems. Some of their applications in energy engineering are briefly described as follows.
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Goldberg (1989, Chapter 4) had applied GAs to natural gas pipeline control since 1983 but his objectives were different from ours. The objective was to minimize the power consumed subject to maximum and minimum pressure and pressure ratio constraints. Square pressure difference across each compressor station was chosen as the control variable. The constraints were incorporated into the system using a quadratic, external penalty method. The study showed that near optimal results were found quickly after searching an infinitesimal portion of the space. Dasgupta (1997) presented an optimal real-time scheduling system for power generation. The objective of scheduling was to generate power to meet the load with minimum cost of generation. The cost was the sum of start-up, banking, and expected running costs subject to the demand and spinning reverse. The constraints included minimum down time and up time of generating units. A penalty function was added to the objective function to penalize a solution when a constraint is violated. The study suggested that the GA based method is a feasible alternative approach for this type of optimization problem. However, a drawback of this approach is the computational time needed to evaluate the population in each generation. Lu et al. (2005) presented a GA approach to optimize the in-building section of the centralized heating, ventilation and air-conditioning (HVAC) systems that consist of indoor and chilled water loops. The approach aimed to optimize energy consumption of the overall system rather than those of individual components, and unknown parameters were calculated using mathematical models and a neuro-fuzzy inference system. The objective was to minimize the sum of power consumption for chillers, pumps and fans, which was subject to five constraints. Two of these constraints involve structures that make it difficult to obtain their derivatives, rendering it difficult to solve them with traditional techniques. The GAs handled these constraints with penalty functions; and the results showed that the method significantly improved system performance. Osman et al. (2004) proposed a methodology to optimize power flow as a step towards optimal operation and control of power generation. The methodology consists of two GAs and handles constraints by introducing repair functions for infeasible solutions. The first GA searches for an initial reference point that violates fewer constraints than a desired precision. The second one uses the feasible reference point to repair the infeasible solution while searching for an optimal solution. The methodology was applied successfully on a six-bus sample power system. Baran et al. (2005) analyzed six different multiobjective evolutionary algorithms in solving an optimal pump-scheduling problem. The six algorithms are: (1) Strength Pareto Evolutionary Algorithm (SPEA), (2)
117
Non-Dominated Sorting GA (NSGA), (3) Non-Dominated Sorting GA 2 (NSGA2), (4) Controlled Elitist Non-Dominated Sorting GA (CNSGA), (5) Niched Pareto GA (NPGA), and Multiple Objective GA (MOGA). The authors recommended SPEA for the problem. The study differed from Goldberg (1989), Dasgupta (1997), Lu et al. (2005) and Osman et al. (2004) in that it handled multiple objectives separately instead of bundling them into a single objective function. Instead of generating a single solution, the optimization process produced a set of solutions that traded off different objectives including minimizing electric cost, maintenance cost, maximum power peak, and level variation in a reservoir. A repairing algorithm was used to ensure that the solutions generated fall within a feasible range.
3. Neural networks for BHP requirement prediction In order to determine an optimal compressor schedule, one of the factors to consider is customer demand which is difficult to predict. Instead, a heuristic method for calculating the required BHP, which can be used as a proxy for customer demand, was adopted in this study. According to Sun et al., 2000 and Uraikul et al., 2000, the BHP requirement for satisfying customer demand can be estimated using flow rates from the two gas stations with the following heuristic equation: BHP requirement ¼ max ð0; 2:77411 ðSt:Louis Flow þ Melfort FlowÞ 1132Þ.
ð1Þ
Since the independent variables of St. Louis and Melfort flow rates in Eq. (1) are unknown at the beginning of a shift, the neural network technique was adopted to predict flow rates of the current shift based on those of the last shift. Fig. 4 shows the overall structure of the system, which includes two main parts: the neural networks and the GAs. The neural networks take past flow rates from Melfort and St. Louis stations and predict the future flow rates. The outputs of the neural networks are then used for calculating the BHP requirements, which represent the customer demand and in turn serve as part of the input for the GAs. The GAs optimize compressor schedule based on several parameters including BHP requirements, compressor technical parameters, penalty factors and GA search parameters. The compressor technical parameters include fuel cost, start-up cost, and horsepower rating. Penalty factors decide the importance of each type of penalty, including wear and tear penalty, oversupply penalty, and noncovered penalty. These penalties will be discussed further in Section 4.1.2. The GA parameters govern the search process.
ARTICLE IN PRESS 118
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 4. Overall structure of the compressor scheduling system.
Neural network models were trained on historical time-series data of flow-rate recorded at the two stations of St. Louis and Melfort. The data were collected in two periods. The first period was from the cold season of December 3rd, 2001 to March 13th, 2002, and the second period was from the hot season of May 28th to August 8th, 2002. Three alternative neural network structures were considered. In all three alternatives, the Melfort and St. Louis models were developed separately. There are eight input variables in each model, each of which represents the flow rate during the period of an hour in the previous shift. However, the alternative neural network structures differ in the number of output variables and number of neural networks involved. The first structure consists of eight output variables, each of which is a flow rate in the current shift. While this model has the advantage that all the flow rates can be calculated in one step, a weakness is that the connection between the input layer and the hidden layer must account for all eight outputs and hence its underlying structure is more complex. The second structure aims at reducing this complexity by using eight separate neural networks; each neural network predicts the flow rate for one of the 8 h in a shift. The disadvantage of this approach is the large number of neural networks involved; altogether, 32 analytical and 16 final models must be built. In the third alternative structure, only the model that predicts 1 h ahead is necessary and this model will be used recursively to predict 8 h ahead. This type of model has low complexity and can adequately explain training data. However, since the models are trained for short term and used for long term forecast, the generalization may be unreliable. In selecting among the three alternative structures, the second alternative was rejected because it involved implementing and managing a large number of neural
Table 1 Models’ validation root mean square errors Station
Training season
Validation season
Number of outputs
RMSE
Melfort Melfort Melfort Melfort St. Louis St. Louis St. Louis St. Louis
Cold Cold Hot Hot Cold Cold Hot Hot
Hot Hot Cold Cold Hot Hot Cold Cold
1 8 1 8 1 8 1 8
113.55 55.93 91.43 74.05 101.53 118.35 225.10 83.92
networks. The remaining two alternatives were implemented and analyzed. Since characteristics of the hot vs. cold seasons are significantly different, the data sets for one season were used to develop models while the other sets were used to validate the models’ generalization. As can be seen in Table 1, the first alternative (highlighted rows) generated lower root mean square errors (RMSE) on the validation data sets for 8 h ahead, and hence it was used for developing the final models. The data sets for both seasons were then combined to build the final models which were incorporated into the GA system. The training results are summarized in Table 2. During training, the data was re-sampled in different ways into two portions for training and testing; and the testing portion ensured not only that the model fitted well to the data but also that it could be generalized to cover unseen data. Both training and testing errors and correlative coefficients were recorded, and plots of model-generated and observed data points for 1 hand 8 h ahead are displayed in Figs. 5 and 6, respectively. As expected, prediction accuracy decreased as the terms of prediction increased. For example, in Table 2, the testing RMSE of the St. Louis model
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
119
Table 2 Training statistic of final models Hours ahead
Melfort
St. Louis
Training
8 7 6 5 4 3 2 1
Testing
Training
Testing
RMSE
CORR_COEFF
RMSE
CORR_COEFF
RMSE
CORR_COEF
RMSE
CORR_COEF
50.88 50.84 49.01 46.87 54.16 42.61 38.43 28.94
0.8392 0.8442 0.8545 0.8648 0.8778 0.8909 0.9115 0.9501
64.79 61.51 59.63 57.14 53.07 48.75 43.31 34.06
0.7485 0.7742 0.7875 0.8033 0.8333 0.8603 0.8906 0.9323
62.97 58.81 57.87 56.43 52.26 47.75 41.52 31.43
0.9152 0.9292 0.9315 0.9343 0.943 0.9526 0.9642 0.9796
70.02 64.25 62.49 60.39 57.22 52.35 45.87 34.06
0.9211 0.9241 0.9284 0.933 0.9398 0.9497 0.9615 0.9785
increased from 34.06 to 70.02 m3/h when the term of prediction increased from 1 to 8 h ahead. This can also be observed from Figs. 5 and 6. The dots in Fig. 5 gather closely to the line while the dots in Fig. 6 are more scattered. The training and testing correlative coefficients (CORR_COEFF in Table 2) of the St. Louis model was higher than those of the Melfort models, which means that the St. Louis model both explains the training data better and has higher generalization ability than the Melfort model. On the other hand, the lowest testing correlative coefficient of the Melfort model was 0.7485, which indicates that the trained model can explain at least 74.85% of the unseen data. This is an acceptable result. The neural network models were implemented in Neur-Online (NOL) Studio (trademark of Gensym Corp., USA) and were imported to the main program as ActiveX components. From the neural network modeling, the flow rates and hence the BHP requirements are obtained, which becomes input to the second subproblem of compressor scheduling. Next, we present representation of the problem of compressor scheduling for GA modeling.
4. Genetic algorithms for optimization of gas compressor schedule 4.1. Problem representation 4.1.1. Chromosome The first step in developing a GA is to encode the problem space of operational schedules. The entire operational schedule was encoded in a single chromosome. There are five controllable compressors operating in a time span of 8 h. The change of compressor states is assumed to occur only at the beginning of each hour. Hence, there are a total of 5 times 8 or 40 bits in a chromosome string. Each bit will have a value of 1 if the compressor is turned on and 0 if the compressor is
turned off. Since there are two choices of value for the first bit, and for each choice of the first bit there are two choices of the second bit, and so on, the total number of possible solutions to this optimization problem is 24041012. It could be difficult for a traditional method to explore such a large search space; hence, GAs, being an efficient heuristic search technique, are employed to analyze the problem. Two alternatives organizations were considered for the order of the chromosome. A chromosome can either contain eight sub-strings for eight periods of five compressors’ states, or five sub-strings of compressors’ states in eight periods. In the first alternative, if a chromosome is represented with an array of bits x[1..40], x[1] to x[5] will represent the states of the five compressors in the first hour, x[6] to x[10] will represent the states of the same compressors in the second hour, and so on. In each sub-string, the first three compressors run on natural gas and the last two run on electricity. In the second alternative, x[1] to x[8] will represent the states of the first compressor in the eights periods of an hour each, x[9] to x[16] will represent the states of the second compressor in the eight periods, and so on. The order difference does not affect fitness of a chromosome; however, it could make a difference in the search since it is more difficult to break up genes that are close together in a crossover process. Since related genes should be placed close to one another, the relationships among genes are examined. Referring to the cost and penalty functions in the next section, a sequence of states of a compressor in a shift influences start-up cost and wear-and-tear penalty. However, calculation of start-up cost relies only on two consecutive periods, and wear-and-tear penalty involves five periods (the calculation of wear-and-tear penalty will be explained in Section 4.1.2). A combination of compressors in a period, on the other hand, influences penalties for oversupply and non-covered demand in that period. Since the non-covered demand penalty is a hard constraint that should not be violated,
ARTICLE IN PRESS 120
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 5. St. Louis Model—predicted vs. observed—1 h ahead.
it is more important than the other costs and penalties. Hence, it is concluded that the relationship between compressors in a period is more important than the sequence of a compressor’s states in eight periods, and the first alternative organization of chromosome order is preferred. The binary coding reduces the number of constraints from seven in the LP representation approach (Sun et al., 2000) to one in the GA approach. The reason for this reduction is explained as follows. In the pipeline problem, compressor states have discrete values of 0 or 1 that can be easily represented as binary bits in a chromosome. However, since discrete parameter values cannot be directly represented in the LP approach, a number of extra inequality constraints were added. Without these extra constraints, the LP system would generate infeasible solutions. For example, in Sun et al.
(2000), to represent the binary state of natural gas compressor i in period p (denoted as Ni, p), a linearization variable xi;p ¼ N i;p N i;p1 was introduced and the following three extra constraints were added. N i;p þ N i;p1 xi;p p1 8i;p ; N i;p N i;p1 þ 2xi;p p0 8i;p ; N i;p þ N i;p1 xi;p X0 8i;p : These three extra constraints are always true in the GA setting and therefore unnecessary. 4.1.2. Fitness evaluation The objective function is specified as the sum of fuel cost, start-up cost, wear and tear penalty (i.e. cost for not letting a compressor run for at least three consecutive hours), oversupply penalty (i.e. cost for
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
121
Fig. 6. St. Louis Model—predicted vs. observed—8 h ahead.
oversupply of gas), and non-covered penalty (i.e. cost for not satisfying customer demand). Since individuals with higher fitness are preferred, the fitness function is expressed as the negation of the objective function.
s(i) d(p)
fitnessðxÞ ¼ ½Fuel_CostðxÞ þ Startup_CostðxÞ
c1(i)
þ Wear_Tear_PenaltyðxÞ þ Oversupply_PenaltyðxÞ þ Non_Covered_PenaltyðxÞ. Let us denote the problem parameters as follows: x[1..40] Array representing a chromosome, where each element x[j], jA[1, 40] represents the state of the (j mod 5)th compressor in the (floor(j/5))th period f(i) One-hour fuel cost for compressor i, iA[0, 4]
h(i)
c2 c3
Start-up cost for compressor i, iA[0, 4] BHP requirement in period p (Proxy for customer demand), pA[1, 8] Horsepower rating of the ith compressor, iA[0, 4] Wear and tear penalty factor for compressor i, iA[0, 4] Oversupply penalty factor Non-covered penalty factor
The partial costs and penalties are calculated as follows:
Fuel cost: The fuel cost is calculated as total fuel cost consumed by all compressors of the two gas stations during a shift. The fuel consumption rates of
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
122
compressors are fixed and different. The cost for running a natural gas compressor is higher than that for an electrical compressor. Fuel_CostðxÞ ¼
40 X
The pseudo-code can be translated into the following equivalent equation: Wear_Tear_PenaltyðxÞ ¼
f ðj mod 5Þx½j.
Start-up cost: The start-up cost is accumulated each time a compressor is turned on from its idle state. Start-up cost for an electrical compressor is higher than that for the natural gas compressors. Start_CostðxÞ ¼
40 X
sðj mod 5Þ maxðx½j x½j 5; 0Þ.
i¼6
Wear and tear: The combination of compressors should be selected in such a way that the number of switching is minimized. After a compressor is turned on, it should be kept in operation for as long as possible. Similarly, when a compressor is shut off, it should be kept inactive for as long as possible. In this study, the minimum operative and inoperative period is set to be 3 h.
c1 ½j mod 5
j¼6
x½j x½j þ 5 þ x½j x½j þ 5 ! 0 mþ2 Y X 2 ð x½j þ 5 kÞ mod 3 .
j¼1
35 X
m¼2 k¼m
This equation is non-linear and cannot be solved using LP. Sun et al. (2000) provided a linear equation to describe the wear and tear factor, which punishes any three consecutive sequences that have a sum of less than three. However, the equation wrongly punishes the free flow solution in which no compressors are running. Another flaw in the approach reported in Sun et al. (2000) is that the wear and tear penalty did not take past states of the compressors into consideration.
Oversupply penalty: This penalty applies when there is surplus gas delivered to the customers. The oversupply amount of gas will be wasted and is assigned a penalty. OverSupply_PenaltyðxÞ
The wear and tear penalty is implemented as follows. The penalty for each compressor is calculated by examining the center of a moving window of five successive periods w[1..5]. If the center w[3] is the same as either of its two predecessors w[2] and w[1], successors w[4] and w[5], or immediate neighbors w[2] and w[4], it should not be punished. Otherwise, the penalty will be calculated based on the center’s immediate neighbors. If we do not have enough information to determine whether or not we should penalize a period, we will not punish the solution. For example, if information on the previous shift is not available then the first period of the current shift should never get penalized for wear and tear. Let us denote st[1..8] as the array of the states of a compressor during a shift. In the current implementation, the eight periods 1–8 are extended to 10 periods 0–9 where states st[0] ¼ st[1] and st[9] ¼ st[8]. The pseudo-code for calculating the wear and tear penalty is as follows, where p refers to index in [1..8].
For period p from 1 to 7 If (not ((st[p] ¼ st[p1] and st[p] ¼ st[p2]) or (st[p] ¼ st[p1] and st[p] ¼ st[p+1]) or (st[p] ¼ st[p+1] and st[p] ¼ st[p+2]))) then penalty ¼ penalty+c1*(|st[p]st[p1]| +|st[p]st[p+1]|) End if End for
¼ c2 max 0;
40 X
!!2 hðj mod 5Þx½j dðpÞ
.
j¼1
Penalty for non-covered demand: This penalty applies when the gas supply cannot satisfy customer demand. If the supply is sufficient, no penalty is applicable; otherwise, the part of demand not covered by supply will be punished. The penalty factor c2 of this penalty should be much higher than the penalty factor c1 of the oversupply penalty because if demand is not covered, customer satisfaction cannot be ensured. Non_Covered_PenaltyðxÞ ¼ c2 max 0; dðpÞ
40 X
!!2 hðj mod 5Þx½j
.
j¼1
The constraints in the GA representation are more compact than those in the LP approach. In the GA approach, all constraints are bundled in one value while there are eight inequalities for the eight periods in the LP approach. 4.1.3. Genetic algorithm process The simple GA that was introduced briefly in Section 2.3 was used to solve the compressor-scheduling problem. The algorithm started with a small number of random individuals in the initial population. At each
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
generation, the better individuals were chosen with higher probabilities for creating successors. The best individuals from each generation were kept. The runtime of the GA for 1000 generations took only a few minutes. This efficiency in runtime can possibly be explained as follows: (1) the GA approach only examines a small number of individuals each time, (2) it only searches in some parts of the search space that are deemed to be closer to the global optima instead of searching the entire problem space as in exhaustive search, and (3) the time for calculation and comparison of fitness function is likely to be less than the time needed to solve equations.
a value of 30 or 50, the number of generations varied from 100 to 1000, the crossover rate had values of 0.6, 0.8 or 1.0 and the mutation rate had values of 0.0005, 0.004, and 0.03. The GAs were executed from three to five times with each combination of these values. Although the results are different in each run due to randomness in the initial population, the following observations can be made.
4.2. Discussions This section discusses the GA search process from two pespectives. First, the impact of parameter values on the search process of the GA will be examined. Secondly, the penalty factor will be discussed. The GAs were executed with the fixed parameters of fuel cost, start-up cost, and horsepower rating. The values for the parameters are listed in Table 3. Other parameters may vary. As can be seen from Table 3, the two electrical compressors (E1 and E2) have the same technical specification. Similarly, so do the two natural gas compressors (N2 and N3) at Melfort station. The natural gas compressors have higher horsepower rating but also require more fuel to operate. 4.2.1. Genetic algorithm parameters In order to examine the impact of parameter values on the search process, the parameter values in the GA were varied as follows. The population size was assigned
Table 3 Fixed parameters Compressor
Fuel cost
Start-up cost
Horsepower rating
N1 N2 N3 E1 E2
400 350 350 100 100
120 120 120 200 200
650 600 600 250 250
123
The change in population size and crossover rate did not affect the genetic search in any way. When the BHP requirements were small, the search converged quite quickly. However, when the BHP requirements increased, the best solution usually appeared later. Hence, a large maximum number of generations are always recommended. The mutation process is needed to recover useful genetic material that was accidentally lost over generations. Hence, as can be observed in Fig. 7, when the mutation rate was too small (0.0005), the search was trapped and could not return to its initial good solutions in later generations. The initial best solution has the fitness values in the range of 80–90 but since the fitness drop to 50 and under in a bad transition, it cannot regain its initial values. On the other hand, if the mutation rate is too high (0.03) as in Fig. 8, the search oscillated randomly. Finally, when the rate is right, the fitness improved through generations (see Fig. 9). The charts on the left of Figs. 7–9 represent the scaled fitness of the best solutions through generations and the charts on the right represent the average solutions.
4.2.2. Penalty factors The penalty factors or coefficients have significant impact on the quality of solutions. Among the factors, the penalty factor for non-covered demand is the most important, because it ensures customer satisfaction and represents a hard constraint that should not be violated in the final best solution. While some infeasible solutions in the population can provide diversity, the number of infeasible solutions should not be too large.
Fig. 7. Effect of a small mutation rate.
ARTICLE IN PRESS 124
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
Fig. 8. Effect of a large mutation rate.
Fig. 9. Effect of a good mutation rate.
It is observed that the penalty factor for non-covered demand (c3) should be proportional to the BHP requirements. If the BHP requirements are small, the feasible search space is large, so even if the penalty is small, there are still some good feasible solutions in the population. However, when the BHP requirements are large, the penalty should be stricter in order to maintain feasibility of population. Our experiment showed that when the BHP requirements are 2000 or above, the penalty factor for noncovered demand must be 3,000,000 or above in order to have a feasible solution. However, even with such a high factor, the number of feasible solutions in all generations was low. Nevertheless, the average number of constraint violations of each chromosome seemed to reduce from 4 to 6 each at early generations to 1–3 in later generations. The disadvantage of a high non-covered penalty factor is that it will overshadow the importance of the wear and tear and other cost and penalties. In the case where there are mixed low and high BHP requirements (e.g. low in periods 1–4 and high in periods 5–8), this problem can produce costly solutions. For example, let us look at a solution when the BHP requirement is high as shown in Table 4. The highest BHP requirement in the eight periods is 2000, which means the non-covered penalty factor is high. This overshadows for example, the wear and tear cost. Since E1 and E2 are of the same capacity, the solution is not optimal because of the unnecessary back and forth switches between E1 and E2 from period 3 to 4 and from period 5 to 6. The solution
Table 4 A genetic algorithm solution Period
BHP requirement
N1
1 2 3 4 5 6 7 8
800 1000 1100 1500 1600 2000 1300 1100
On
On On On On On
N2
On On On On On
N3
E1
E2
On On On
On
On On On
On On On On
On On On On
On On
is better if E1 remains off and E2 remains on in periods 4 and 5. Hence, each period should have its own penalty factor for non-covered demand instead of sharing the same penalty factor as in the current implemented system. This can be considered in future work.
5. Results: comparison with past solutions In this section, the solutions from the three methods of expert system, LP, and GAs will be comprared in two different cases with different sets of BHP requirements. As all the solutions satisfy the hard constraint of ensuring customer satisfaction, the comparison was made in terms of wasted energy and compressor switches between consecutive periods. As mentioned in Section 2.1, the evaluation version of the HYPER LINDO software package for LP could
ARTICLE IN PRESS H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
handle a maximum of only 200 constraints, and only six periods can be optimized at each run. To make a reasonable comparison, the BHP requirements from Sun et al. (2000) were used for periods 1–6 and two random BHP requirements were used for periods 7 and 8. Table 5 displays the solutions from the three methods of expert system, LP, and GAs. The GA solution wasted less energy than the expert system but is not as efficient as the LP solution. For example, in the fifth period, the expert system used all five compressors while the GA turned off an electrical compressor. The LP approach went a step further when it substituted a gas compressor in the GA solution with a less powerful electrical compressor. In the first six periods, the total horsepower supplied with the schedules generated by expert system, LP and GA are 10950, 8700, and 9750, respectively, while the necessary horsepower was only 8000. Given the discussion in Section 4.2.2, the wasted resources can be accounted for as follows: since the BHP requirement is high, the non-covered demand penalty factor was too high and hence, other costs and penalties were neglected. However, in terms of wear and tear penalty, the expert system provided the most efficient solution because in the first six periods, the BHP requirements increase monotonously and hence compressors were only added and not turned off. The only penalty that the expert system was subjected to is when compressor N3
125
was turned off in period 7 after only two running periods. The GA solution was inefficient in terms of wear and tear cost, because compressor switches were frequent in almost every period of operation. On the other hand, the LP seemed to ignore the wear and tear factor as they allowed the OR operator in its set of solutions. As an example, if the combinations of the first two periods are {E1, N2} and {E1, E2, N2}, then there are no switches; but if they are {E1, N2} and {E1, E2, N3}, then there is a switch between two equivalent natural gas compressors N2 and N3. Another experiment was conducted; this time the BHP requirements were predicted based on previous flow using neural networks. No results from LP are available for this set of BHP requirements. The results from GA are shown in Table 6. Again in this case, we observe that the solution from GA wasted less resource than that from the expert system. The expert system experienced a sudden shift in the amount of energy generated when the BHP requirement crossed the 800 level. In periods 1 and 2, when the BHP requirements were close but still less than the upper bound of 800 HP, the expert system was very efficient. However, it was least effective when it barely crossed the boundary of 800 HP because the total energy generated by the expert system approach was 1500 HP, which means (1500–800) or 700 units of horsepower was
Table 5 Comparison of the three approaches Period
BHP requirement
Expert system
Linear programming
Genetic algorithm
1 2 3 4 5 6 7 8
800 1000 1100 1500 1600 2000 1300 1100
E2, E2, E2, E2, E2, E2, E2, E2,
(E1|E2), (N2|N3) E1, E2, (N2|N3) E1, E2, (N1|N2|N3) E1, E2, N1, (N2|N3) E1, E2, N1, (N2|N3) (E1|E2), N1, N2, N3
N1, N3, N2, N1, N1, N1, N1, N1,
N1, N1, N1, N1, N1, N1, N1, N1,
N2 N2 N2 N2, N2, N2, N2, N2
E1 N3, E1 N3, E1 E1
N3, E2 E1, E2 N3, E2 N2, E1 N2, N3, E1 N2, N3, E2 N3, E1, E2 N2, N3, E1
Table 6 Genetic algorithm solution in the case of medium BHP requirement Period BHP requirements Expert system solution Genetic algorithm solution Total energy generated by ES Total energy generated by GA 1 2 3 4 5 6 7 8 Total
796.09 799.73 803.98 801.71 811.94 813.71 818.88 826.88 6472.92
E2, E2, E2, E2, E2, E2, E2, E2,
N1 N1 N1, N1, N1, N1, N1, N1,
N2 N2 N2 N2 N2 N2
N1, N1, N2, N1, N3, N2, N2, N3,
E1, E2 N3 E1 N3, E1 E2 N3 E1 E1, E2
900 900 1500 1500 1500 1500 1500 1500 10800
1150 1250 850 1500 900 1200 850 1100 8800
ARTICLE IN PRESS 126
H.H. Nguyen, C.W. Chan / Engineering Applications of Artificial Intelligence 19 (2006) 113–126
wasted. In terms of wear and tear cost, the GA solution was also inefficient in this experiment because of frequent switching.
population replacement or De Jong’s crowding model (Goldberg, 1989).
Acknowledgment 6. Conclusions and future work The genetic algorithms (GAs) in this study provided an efficient way to search for a compressor operation schedule. The running time for 1000 generations of 30 individuals each is less than 5 min which is significantly shorter than that of the linear programming approach. However, it is difficult to determine if a near optimal solution has been found or whether the search has been trapped in a local minimum. The quality of the solution also depends heavily on the initial random population. The fact that GAs use payoff information, not derivatives or other auxiliary knowledge, gives it more flexibility than other traditional methods. GAs can handle non-linearity in the same way as linearity. However, this heuristic also involves a weakness in that some information of the function to be optimized is lost in the process of calculating payoffs. The approach to handle constraints in this study is not yet efficient. It depends largely on a number of penalty coefficients that need to be tuned. In future research, the use of repairing technique to repair infeasible solutions will be considered. A simple repairing algorithm could be as follows. Until the demand is met, keep turning on more compressors, which are selected at random. In the case of small feasible search space such as when high BHP level is required, the repairing approach is especially helpful since it ensures that a feasible solution is found. Another possible direction for future work is to add other constraints such as maximum operating hours or real-time failure. The genetic algorithms can also be extended with more complex features such as elitist
The authors would like to acknowledge the generous support of a Strategic Grant from the Natural Sciences and Engineering Research Council of Canada. References Baran, B., von Lucken, C., Sotelo, A., 2005. Multi-objective pump scheduling optimization using evolutionary strategies. Advances in Engineering Software 36, 39–47. Coello, C.A., 2002. Theoretical and numerical constraint-handling techniques used with evolutionary algorithms: a survey of the state of the art. Computer Methods in Applied Mechanics and Engineering 191 (11–12), 1245–1287. Dasgupta, D., 1997. Optimal Scheduling of thermal power generation using evolutionary algorithms. In: Dasgupta, D., Michalewicz, Z. (Eds.), Evolutionary Algorithms in Engineering Applications. Springer, Berlin, Heidelberg, pp. 317–328. Dasgupta, D., Michalewicz, Z., 1997. Evolutionary algorithms—an overview. In: Dasgupta, D., Michalewicz, Z. (Eds.), Evolutionary Algorithms in Engineering Applications. Springer, Berlin, Heidelberg, pp. 3–28. Goldberg, D.E., 1989. Genetic Algorithms in Search, Optimization, and Machine Learning. Addition-Wesley Publishing Company Inc., Reading, MA. Lu, L., Cai, W., Xie, L., Li, S., Soh, Y.C., 2005. HVAC System optimization—in-building section. Energy and Buildings 37, 11–22. Osman, M.S., Abo-Sinna, M.A., Mousa, A.A., 2004. A solution to the optimal power flow using genetic algorithm. Applied Mathematics and Computation 155, 391–405. Sun, C.K., Uraikul, V., Chan, C.W., Tontiwachwuthikul, P., 2000. An integrated expert system/operation research approach for the optimization of natural gas pipeline operations. Engineering Applications of Artificial Intelligence 13, 465–475. Uraikul, V., Chan, C.W., Tontiwachwuthikul, P., 2000. Development of an expert system for optimizing natural gas pipeline operations. Expert Systems with Applications 18, 271–282.
