EVOLUTIONARY COMPUTATION AND
CONSTRAINED OPTIMIZATION
Ruhul Sarker, Thomas P. Runarsson* and Charles Newton
School of Computer Science
UC, University of New South Wales, ADFA
Northcott Drive, Canberra, ACT 2600, Australia
<ruhul@cs.adfa.edu.au>
*Current Address: Dept. of Mechanical Engineering, University of Iceland, Iceland
Abstract
We consider a class of constrained nonlinear integer programs, which arise in
manufacturing batch sizing problems. We investigate the use of genetic algorithms
(GAs) for solving these models. In this paper, a new penalty function based GA
method is proposed. Both binary and real coded genetic algorithms with six other
penalty functions are developed and compared with the proposed method. The real
coded genetic algorithm works well for all penalty functions compared to binary
coding and the new method shows the best performance. Numerical examples are
provided and computational experiences are discussed.
1. Introduction
Consider a constrained nonlinear programming problem.
Problem P1
Minimize f(x)
Subject to g(x) 0
h(x) = 0
xX
where g is a vector function with components g1, ……, gm and h is a vector function with
components h1, ……, hl. Here f, g1, ……, gm, h1, ……, hl are functions on En and X is a
nonempty set in En. The set X represents simple constraints that could be easily handled
explicitly, such as lower and upper bounds on the variables.
In contrast to the unifying role of the simplex method in linear programming, there
exists no unique approach to nonlinear optimization. Numerous algorithms have been
proposed to solve the nonlinear programming problem such as method of steepest decent [1]
for solving unconstrained optimization, method of Lagrange multipliers [2] for solving
equality constrained optimization, and the gradient projection method [3]. The gradient
projection method provides a good procedure for moving along the boundary of linearly
constrained regions. In order to solve a general nonlinear constrained programming
problem, many different techniques have been proposed. Courant [4] initially suggested the
penalty function, which studies constrained problems via related unconstrained ones. In the
penalty function method, a penalty term is added to the objective function for any violation
of the constraints.
1
Similar to the penalty function, Carroll [5] proposed the barrier function approach in
which a barrier term that prevents the points generated from leaving the feasible region is
added to the objective function. This method can be used for inequality constrained
problems. Abadie and Carpenter [6] generalized the reduced gradient methods to nonlinear
constraints. The method is an extension of the Wolfe [7] algorithm to accommodate both
nonlinear objective functions and nonlinear constraints. These methods have been studied in
detail in the past and have been found to have weaknesses. In the penalty and barrier
function methods, the unconstrained subproblem becomes extremely ill-conditioned for
extreme values of the penalty/barrier parameters. Reduced gradient methods have
difficulties following the boundary of high nonlinear constraints. These methods rely on
local gradient information; the optimum is the best in the neighborhood of the current point.
To enhance the reliability and robustness of the search technique in constrained
optimization, a GA based search method incorporating several penalty functions is tested in
this paper.
We consider a manufacturing batch sizing problem. A single raw material
unconstrained manufacturing batch sizing model was developed by Sarker and Khan [8].
Later, this model was modified by Sarker and Newton [9] to incorporate the transportation
module for finished products delivery. In this paper, we considered a problem similar to
Sarker and Khan [8] but with multiple raw materials and several constraints to make the
situation more realistic. The new problem forms a constrained nonlinear integer program.
Considering the complexity of solving such a model, we investigated the use of genetic
algorithms (GAs). We used both binary coding and real coded genetic algorithms with
different crossovers and mutations. We also proposed a new penalty function based GA
method. This method provides better solutions than the existing methods. It is found in the
experimentation that the real coded genetic algorithm works well for this batch sizing
problem. The detailed computational experiences are presented.
The organization of the paper is as follows: Following this introduction, this paper
presents a brief introduction on penalty function methods and then on genetic algorithms. In
section 4, the batch sizing problem is described, and its mathematical formulation is presented.
In section 5, the different transformation methods are provided. The model solving approach
using GAs is implemented in section 6. The computational experiences are presented in section
7. Finally, conclusions are provided in section 8.
2. Penalty Function Method
The penalty function method has been well known for solving constrained optimization
problems for decades. The approach converts the problem into an equivalent unconstrained
problem and then solves it using a suitable search algorithm. Two basic types of penalty
functions exist: exterior penalty functions, which penalize infeasible solutions, and interior
penalty functions, which penalize feasible solutions. We will discuss only the exterior
penalty functions in this paper as the implementation of interior penalty functions is
considerably more complex for multiple constraint cases [10].
Three degrees of exterior penalty functions exist: (i) barrier methods in which no
infeasible solution is considered, (ii) partial penalty functions in which a penalty is applied
near the feasibility boundary, and (iii) global penalty functions that are applied throughout
the infeasible region [11]. In general a penalty function approach places the constraints into
the objective function via a penalty parameter in such a way that it penalizes any violation
of the constraints. A general form of the unconstrained problem is as follows:
Problem P2
Minimize fp(x) = f(x) + (x)
2
Subject to x X
where > 0 is a large number, and (x) is the penalty function.
We define the problem P2 as a penalty problem and the objective function of P2 as the
penalized objective function (fp(x)) which is the simple sum of the unpenalized objective
function and a penalty (for a minimization problem). A suitable penalty function incurs a
positive penalty for infeasible points and no penalty for feasible points. The penalty function
is usually of the form:
(x) =
m
[maximum{0, gi(x)}]p +
i 1
l
| hi(x)|p
…..
…..
……… (1)
i 1
where p is a positive integer.
If the penalty parameter = 0, the problem P2 forms an unconstrained optimization
problem. Normally an independent penalty parameter is imposed for the violation of each
constraint i, instead of having a constant for all constraints.
The solution to the penalty problem can be made arbitraily close to the optimal solution
of the original problem by choosing sufficiently large. However, if we choose a very large
and attempt to solve the penalty problem, we may get into some computational difficulties
of ill-conditioning [12]. With a large , more emphasis is placed on feasibility, and most
procedures will move quickly toward a feasible point. Even though this point may be far
from optimal, premature termination could occur.
As a result of the above difficulties associated with large penalty parameters, most
algorithms use penalty functions that employ a sequence of increasing penalty parameters
[k+1 = . With each new value of the penalty parameter, an optimization technique is
employed, starting with the optimal solution corresponding to the previously chosen
parameter value.
It is very difficult to choose an appropriate [13]. A lower value of means a higher
number of iterations is required to solve the penalty problem. A number of methods for
selecting are discussed in a later section.
3. Introduction to GA
During the last two decades there has been a growing interest in algorithms which are
based on a principle of evolution- survival of the fittest. A common term, accepted recently,
refers to such techniques as evolutionary computation (EC) methods. The best known
algorithms in this class include genetic algorithms, evolutionary programming, evolution
strategies, and genetic programming. There are also many hybrid systems which incorporate
various features of the above paradigms, and consequently are hard to classify; they are just
referred to as EC methods [14].
The methods of EC are stochastic algorithms whose search methods model some
natural phenomena: genetic inheritance and darwinian strife for survival. GAs follow a stepby-step procedure that mimic the process of natural evolution, following the principles of
natural selection and "survival of the fittest". In these algorithms a population of individuals
(potential solutions) undergoes a sequence of unary (mutation type) and higher order
(crossover type) transformations. These individuals strive for survival. A selection scheme,
biased towards fitter individuals, selects the next generation. This new generation contains a
higher proportion of the characteristics possessed by the "good" members of the previous
generation, in this way good characteristics are spread over the population and mixed with
3
other good characteristics. After some number of generations, the program converges and
the best individuals represent a near-optimum solution. The GA procedure is shown on the
next page.
begin
t0
initialize Population (t)
evaluate Population (t)
while (not terminate-condition) do
begin
tt+1
select Population (t) from Population (t-1)
alter Population (t)
evaluate Population (t)
end
end
4. Batch Sizing Problem (BSP)
We consider a batch manufacturing environment that processes raw materials procured
from outside suppliers to convert them into finished products for retailers. The
manufacturing batch size is dependent upon the retailer’s sales volume (/market demand),
unit product cost, set-up cost, and inventory holding cost. The raw material purchasing lot
size is dependent upon the raw material requirement in the manufacturing system, unit raw
material cost, ordering cost and inventory holding cost. Therefore, the optimal raw material
purchasing quantity may not be equal to the raw material requirement for an optimal
manufacturing batch size. To operate the manufacturing system optimally, it is necessary to
optimise the activities of both raw material purchasing and production batch sizing
simultaneously, taking all operating parameters into consideration. The detail problem
description can be found in [8 & 9].
The unconstrained models under various situations are presented in Sarker and Khan [8]
and Sarker and Newton [9], and the single raw material constrained model was recently
developed by Sarker and Newton [15]. The mathematical formulation for the multiple raw
materials constrained model is presented in this section. The notations used in developing the
model are as follows:
Dp
Pp
Qp
Hp
Ap
ri
=
=
=
=
=
=
Di
Qi
Ai
Hi
PRi
Q*i
=
=
=
=
=
=
demand rate of a product p, units per year
production rate, units per year (here, Pp>Dp)
production lot size
annual inventory holding cost, $/unit/year
setup cost for a product p ($/setup)
amount/quantity of raw material i required in producing one unit
of a product
demand of raw material i for the product p in a year, Di = riDp
ordering quantity of raw material i
ordering cost of a raw material i
annual inventory holding cost for raw material i
price of raw material i
optimum ordering quantity of raw material i
4
xp
=
L
T
mp
ni
=
=
=
=
si
TC
=
=
shipment quantity to customer at a regular interval
(units/shipment)
time between successive shipments = x/Dp
cycle time measured in years = Qp/Dp
number of full shipments during the cycle time = T/L
number of production cycles which consume one lot of raw
material i = Qi/Qp
space required by one unit of raw material i
total cost of the system
The mathematical model is presented below:
Dp
A
Minimize TC
( Ap + i )+
mp xp
i ni
mp xp
2
((
xp
Dp
Dp
+ ni - 1) ) + 1) H p + r i H i (
Hp
2
Pp
Pp
i
Subject to:
n r m
p
x p si Raw _ s _ Cap
n r m
p
x p si Min _ Truck _ Load
i i
i
i i
i
m p x p s p Finish _ S Cap
mp and ni are integers and greater than zero.
This is clearly a nonlinear-integer program. The first constraint indicates that the storage
space required by the raw materials must be less than the space specified for raw materials.
The second constraint represents the lower limit truck load for transportation and the third
constraint is for the finished products storage capacity. If we assign i = 1, then the model
represents a single raw material batch sizing problem.
5. Transformation Methods
The transformation method, from a constrained to an unconstrained problem, in
evolutionary algorithms is almost the same as the classical approach but with a slight
difference. The value of k is changed as a function of generation k, but not is equal to k-1.
Some selected transformation techniques used in the evolutionary computation
literature for solving constrained optimization problems are presented below. Most, if not
all, are of the exterior kind which will allow the initial population of solutions to be partially
or completely infeasible. In some techniques, the controlling parameters are updated in a
predetermined manner. The controlling (or penalty) parameters k may be updated each
generation or every n generations. We present our new approach in this section.
5.1 Static Penalties
The method of static penalties [16] assumes that for every constraint we establish a
family of intervals that determine the appropriate penalty coefficient.
For each constraint, create several (n) levels of violation.
5
For each level of violation and for each constraint, create a penalty coefficient ij
(i = 1, 2, …, m+l; j = 1, 2, .., n); higher levels of violation require larger values
of these coefficients.
Start with a random population of individuals (feasible or infeasible).
Evolve the population; evaluate individuals.
It is clear that the results are parameter dependent. It is quite likely that for a given
problem there exists one optimal set of parameters for which the system returns a feasible
near-optimum solution; however, it might be hard to find. A limited set of experiments
reported by Michalewicz [13] indicate that the method can provide good results if violation
levels and penalty coefficients are tuned to the problem.
We use three levels of fixed penalty coefficients for each constraints: 0.5, 0.5x103 and
0.5x106.
5.2 Dynamic Penalties
Joines and Houck [17] proposed dynamic penalties. The authors assumed that k =
(Ck), where C and are constants. The constant p in equation 1 has a value that is greater
than one. A reasonable choice for these parameters is C = 0.5, = p = 2. This method
requires a much smaller number (independent of the number of constraints) of parameters
than the first method. Also, instead of defining several levels of violation, the pressure on
infeasible solutions is increased due to the (Ck) component of the penalty term: towards the
end of the process (for high values of the generation number k), this component assumes
large values.
5.3 Annealing Penalties
The method of annealing penalties, called Genocop II (for Genetic algorithms for
Numerical Optimization of Constrained Problems) is also based on dynamic penalties and
was described by Michalewicz and Attia [18] and Michalewicz [19]. In this system, a fixed
k is used for all constraints of a given generation k, where k = 1/2. An initial temperature
is set and is used to evaluate the individuals. After a certain number of generations, the
temperature is decreased and the best solution found so far serves as a starting point for
the next iteration. This process continues until the temperature reaches freezing point.
Genocop proved its usefulness for linearly constrained optimization problems; it gave a
surprisingly good performance for many functions [19 & 20]. Michalewicz and Attia [18]
proposed the control parameters to be k+1 = 10.k, where 0 = 1 and k remains constant
once it reaches 1x106.
5.4 Adaptive Penalties
Adaptive transformation attempts to use the information from the search to adjust the
control parameters. This is usually done by examining the fitness of feasible and infeasible
members in the current population.
Bean and Hadi-Alouane proposed an adaptive penalty method in 1992 which uses the
feedback from the search process [see 21]. This method allows either an increase or a
decrease of the imposed penalty during evolution as shown below. This involves the
selection of two constants, 1 and 2 (1 > 2 > 1), to adaptively update the penalty function
multiplier, and the evaluation of the feasibility of the best solution over successive intervals
6
of Nf generations. As the search progresses, the penalty function multiplier is updated every
Nf generations based on whether or not the best solution was feasible during that interval.
Specifically, the penalty function is as follows;
k 1
k 1
k / 2
k
if previous N f generations have only infeasible best solutions
if previous N f generations have only feasible best solutions
Otherwise
The parameters were chosen to be similar to the above methods with Nf = 3, 1 = 5, 2 = 10
and 0 = 1. In our experimentation, the values were increased by 10.
5.5 Death Penalty
The death penalty method just rejects infeasible individuals. It can be interpreted as a
truncation selection based on (x) and then followed by a selection procedure based on f(x).
In this method, the initial population must be feasible. In the experiments reported by
Michalewicz [13] it was shown that the method generally gives a poor performance. Also,
the method is not as stable as others; the standard deviation of solutions returned is
relatively high.
5.6 Superiority of Feasible Points
The method of superiority of feasible points was developed by Powell and Skolnick
[22] and is based on a classical penalty approach, with one notable exception. Each
individual is evaluated by the formula:
fp(x) = f(x) + (x) + (t, x)]
where is a constant; however, the component (t, x) is an additional iteration-dependent
function that influences the evaluation of infeasible solutions. The point is that the method
distinguishes between feasible and infeasible individuals by adopting an additional heuristic
rule: For any feasible individual x and any infeasible individual y, fp(x) < fp(y) (i.e., any
feasible solution is better than any infeasible one). The penalties are increased for infeasible
individuals. The method performs reasonably well. However, for some case problems, the
method may have some difficulty in locating a feasible solution [13].
5.7 A New Method
At this point it seems reasonable to re-examine our goals and what transformations
mean for the evolutionary algorithm. In principle these transformations will influence
selection only.
In this paper the commonly used binary tournament selection is applied. In binary
tournament selection two individuals are chosen randomly from the population and the
better of them is allowed to survive to the next generation. Consider this comparison
between two individuals 1 and 2 transformed to take the general form from problem P2:
f1(x) + (x) f2(x) + (x),
for a given penalty parameter > 0. From equation (2) we can write [f1(x) f2(x))/((x) - (x)]. Now suppose c = [f1(x) - f2(x))/((x) - (x)]. This relation holds
the following conditions:
7
1. If1 = 2, common for feasible solutions, then c is infinity. The competition is
decided by objective function.
2. If f2 = f1, the value of c is zero. The competition is decided by penalty function.
3. If f2 < f1 and 2 < 1. c is negative in this case. Comparison is independent of
4. If f2 < f1 and 1 < 2,c is positive with objective function domination.
5. If f2 > f1 and 1 > 2, c is positive with penalty function domination.
Only for the last two conditions is the value for critical to the outcome of the comparison
and this value is c. Evolutionary algorithms are a population based search method.
Population information can be used to estimate the most appropriate value for . This value
would be the one that attempts to balances the number of competitions decided by the
objective and penalty function. This value may be approximated by simulating the binary
tournament selection and computing for each comparison the critical penalty coefficient c.
Then for the actual selection the used will be the average of these values, or more
precisely
1
M
N
i 1, c 0
c
(i ) ,
where N is the number of simulated competitions held and M<N is the number of
competitions for which c > 0. The number of simulated competitions held is at least twice
the size of the population.
6. Solving BSP using Genetic Algorithms
It is generally accepted that any GA to solve a problem must have five basic
components:
problem representation,
a way to create an initial population of solutions,
an evaluation function rating solutions in terms of their "fitness",
genetic operators that alter the genetic composition of parents during reproduction,
and
values for the parameters (population size, probabilities of applying genetic
operators, etc.)
Two types of evolutionary algorithms were tested, for each of the methods discussed in the
previous section: a Simple Genetic Algorithm (SGA) and a Real Coded Genetic Algorithm
(RGA). In the light of the above five components, these two algorithms are discussed briefly
below.
6.1 Problem Representation
The first hurdle to overcome in using GAs is problem representation- we must
represent our problem in such a way that is suitable for handling by a GA, which works with
strings of symbols that in structure resemble chromosomes. The representation often relies
on binary coding. GAs work with a population of competing strings, each of which
represents a potential solution for the problem under investigation. The individual strings
within the population are gradually transformed using biological based operations. For
8
example, two strings might 'mate' and produce an offspring that combines the best features
of its parents. In accordance with the law of the survival of the fittest, the best-performing
individual in the population will eventually dominate the population.
6.2 Initialize the Population
There are no strict rules for determining the population size. Larger populations ensure
greater diversity but require more computing resources. We use a population size of 50, a
commonly used size (Back, 1996, page 123). Once the population size is chosen, then the
initial population must be randomly generated. For binary coding, the random number
generator generates random bits (0 or 1). Each individual in the population is a string of n
bits. For real coding, the actual numbers are generated randomly.
6.3 Calculate fitness
Now that we have a population of potential problem solutions, we need to see how
good they are. Therefore, we calculate a fitness, or performance, for each string. The fitness
function is the penalty objective (fp) function in our case. For binary coding, each string is
decoded into its decimal equivalent. This gives us a candidate value for the solution. This
candidate value is used to calculate the fitness value.
6.4 Selection and Genetic Operators
Genetic algorithms work with a population of competing problem solutions that
gradually evolve over successive generations. Survival of the fittest means that only the
best-performing members of the population will survive in the long run. In the short run we
merely tip the odds in favour of the better performers; we do not eliminate all the poor
performers. This can be done in a variety of ways. We use a binary tournament selection,
where two individuals are selected at random from the population pool and the better one is
allowed to survive to the next generation.
The chromosomes, which survive the selection step undergo genetic operations,
crossover and mutation. Crossover is the step that really powers the GA. It allows the search
to fan out in diverse directions looking for attractive solutions and permits two strings to
'mate'. This may result in offspring that are 'fitter' than their parents. Crossover is
accomplished for binary coding in four small steps as follows.
two potential parents are randomly chosen from the population;
generate the crossover probability to answer 'yes, perform crossover' or 'no';
if the answer is 'yes', randomly choose a cutting point; and
cut each string at the cutting point and create two offspring by gluing parts of the
parents.
Usually, the classical 1-point crossover is used. Sometimes more than one cutting point
is used. For example, in case of two 2-point crossovers, the first segment of the first
chromosome is followed by the second segment of the second chromosome.
We use a Heuristic crossover for real coding GAs. This is a unique crossover for the
following reasons:
it uses values of the penalty objective function in determining the direction of the
search,
9
it produces only one offspring, and
it may produce no offspring at all. There is a 80% probability that a variable gets
“crossed”.
The operator generates a single offspring x3 from two parents x1 and x2 according to the
following rule:
x3 = rounded[r.(x2 - x1) + x2,]
where r is a random number between 0 and 1, and the parent x2 is no worse than x1.
Actually, when the generated value x3 is out of bounds (upper and lower bounds) it is set
equal to the corresponding violated bound. Also, note that x1 and x2 are “parent variables”
not parents. Also x3 is rounded after crossing to integer. This is not required for binary
coding. Rounding is also performed after the mutation operator.
Mutation introduces random deviations into the population. For binary coding, mutation
changes a 0 into a 1 and vice versa. Mutation is usually performed with low probability,
otherwise it would defeat the order building being generated through selection and
crossover. Mutation attempts to bump the population gently into a slightly better course.
Nonuniform mutation is also used in real coding GA, where
x
t 1
k
x kt (t , u (k ) x k )
t
x k (t , x k l (k ))
if a random binary digit is 0
if a random binary digit is 1
for k = 1, …, n. The u(k) and l(k) are upper and lower domain bounds of the variable xk. The
function (t,y) returns a value in the range [0,y] such that the probability of (t,y) being
close to 0 increases as t increases (t is the generation number). This property causes this
operator to search the space uniformly initially (when t is small), and very locally at a later
stage. Michalewicz et al. [20] used the following function for experimentation:
t
(t , y ) y.r.(1 ) b
T
where r is a random number from (0..1), T is the maximal generation number, and b
(b=6.r+1 was used) is a system parameter determining the degree of nonuniformity.
6.5 Parameters of GA
All GA runs have the following standard characteristics:
Probability of crossover: 1.0
Probability of mutation: 1/(string length of the chromosome)
Population size: 50
Number of generations in each run: 200
Number of independent runs: 300
for SGA: binary coding, two point crossover, and bit-wise mutation.
for RGA: heuristic crossover and non-uniform mutations
7. Computational Results
10
The mathematical model presented in Section 4 is solved for a test problem of a single
product with three raw materials using seven different penalty functions based on binary and
real coded GAs described in Section 5. The data used for the model are as follows:
Sample Data: Dp = 4x106 ; Pp= 5x106 ; Ap= $ 50.00 ; Hp= $ 1.20 , xP = 1,000 units, Sp = 1.00.
Other data are as follows:
Ai
Hi
si
ri
Single Raw Material
For i = 1
3,000
1.00
1.00
1.00
Multiple
Raw material 1
3,000
1.00
1.00
0.5
Raw
Raw material 2
2,500
1.20
1.10
0.2
Materials
Raw material 3
4,600
1.50
1.20
0.3
The right hand sides (RHS) for different test problems are given below:
RHS
Raw_S_Cap
Truck_Min_load
Fin_S_Cap
Problem#1
400,000
200,000
20,000
Problem#2
300,000
100,000
18,000
Problem#3
200,000
40,000
17,000
Problem#4
150,000
30,000
15,000
Problem#5
130,000
25,000
13,000
For single raw material with RHS similar to problem#1, the best solution found was: m
= 13 and n = 10 with TC = $1.8483E5. For multiple raw materials case, the best solution
found for problem#1 was: m = 15, n1 = 15, n2 = 19, and n3 = 19 with TC = $3.3471E5 and
for problem#5 was: m = 12, n1 = 8, n2 = 11, and n3 = 12 with TC = $4.2837E5.
A sample detail comparison of results, for single raw material, provided by seven
different penalty functions are shown in Table 1 for real coded GAs and in Table 2 for
binary coded GAs. The minimum, presented in column 2, means the best objective function
value obtained in any of 300 runs and maximum, presented in column 6, indicates the worst
value. The statistics (like mean, standard deviation and median) of the objective function
values based on 300 runs are also produced to sense the variability of the target value. The
right most columns represent the number of times we hit the optimal solution out of 300
runs, expressed in percentage.
The detail sample comparison of results, for multiple raw materials case problem#1,
provided by seven different penalty functions are shown in Table 3 for real coded GAs and
in Table 4 for binary coded GAs.
Table 1: Results of single raw material case using real coded GAs
Method
New/ Proposed
Death Penalty
Superiority of
Feasible points
Dynamic Penalty
Annealing Penalty
Static Penalty
Adaptive Penalty
Minimum
1.8483
1.8483
1.8483
Mean
1.8483
1.8483
1.8492
Std. Dev.
0.0000
0.0000
0.0047
Median
1.8483
1.8483
1.8483
Maximum
1.8483
1.8483
1.8868
% times Opt. Sol.
100
100
96
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
0.0000
0.0000
0.0000
0.0000
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
1.8483
100
100
100
100
Table 2: Results of single raw material case using binary coded GAs
Method
New/ Proposed
Death Penalty
Superiority of
Feasible points
Dynamic Penalty
Annealing Penalty
Minimum
1.8483
1.8483
1.8483
Mean
1.8484
1.8558
1.9266
Std. Dev.
0.0010
0.0320
0.1521
Median
1.8483
1.8483
1.8868
Maximum
1.8651
2.2653
3.5813
% times Opt. Sol.
100
3
22
1.8483
1.8483
1.8484
1.8484
0.0014
0.0014
1.8483
1.8483
1.8651
1.8651
99
99
11
Static Penalty
Adaptive Penalty
1.8483
1.8483
1.8484
1.8484
0.0010
0.0014
1.8483
1.8483
1.8651
1.8651
100
99
Table 3: Results of multiple raw materials case problem#1 using real coded GAs
Method
New/ Proposed
Death Penalty
Superiority of
Feasible points
Dynamic Penalty
Annealing Penalty
Static Penalty
Adaptive Penalty
Minimum
3.3471
3.3471
3.3471
Mean
3.3471
3.3471
3.3472
Std. Dev.
0.0000
0.0000
0.0002
Median
3.3471
3.3471
3.3472
Maximum
3.3471
3.3471
3.3482
% times Opt. Sol.
100
100
35
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
0.0000
0.0000
0.0000
0.0000
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
3.3471
100
100
100
100
Table 4: Results of multiple raw materials case problem#1 using binary coded GAs
Method
New/ Proposed
Death Penalty
Superiority of
Feasible points
Dynamic Penalty
Annealing Penalty
Static Penalty
Adaptive Penalty
Minimum
3.3472
3.3471
3.3561
Mean
3.3773
3.4478
3.4420
Std. Dev.
0.0210
0.0881
0.0535
Median
3.3853
3.4151
3.4314
Maximum
3.4344
3.6287
3.6462
% times Opt. Sol.
0
1
0
3.3471
3.3471
3.3471
3.3471
3.3773
3.3802
3.3789
3.3774
0.0232
0.0256
0.0219
0.0223
3.3795
3.3933
3.3850
3.3961
3.5313
3.6076
3.4661
3.4270
1
1
1
4
The results, in terms of percentage of times the optimal solution was obtained, are
always better with real coded GAs as compared to binary coded GAs for all seven penalty
functions considered in this paper. All the penalty functions, except Superiority of Feasible
Points, worked extremely well with real coded GAs for the above two test problems. The
results are worse with tighter constrained problems which are demonstrated in Table 5. The
Death Penalty method along with Superiority of Feasible Points performed badly for tighter
constraints. However, all seven penalty functions, with real coding GAs, provide unique
optimal solutions for all the test problems. These solutions are acceptable irrespective of
how many times the optimal solutions are generated in a single run. The binary coded GAs
failed to produce any optimal solution for the tighter constrained test problems. The
proposed method seems a bit superior for the tighter constrained problems.
Table 5: Percentage times optimum solution obtained by different methods
Method
New/ Proposed
Death Penalty
Superiority of
Feasible points
Dynamic Penalty
Annealing Penalty
Static Penalty
Adaptive Penalty
Problem#1
100
100
36
Problem#2
100
90
30
Problem#3
80
9
7
Problem#4
64
6
3
Problem#5
49
5
2
100
100
100
100
100
100
100
100
76
74
78
73
52
54
58
58
42
41
37
38
8. Conclusions
We considered a single product- multiple raw materials joint batch sizing problem. The
mathematical programming model of this problem forms a constrained nonlinear integer
program. It is difficult to solve such a complex model using the existing algorithms. In this
12
paper, we investigate the use of genetic algorithms (GAs) for solving such a model with one
product and single or three raw materials. We used both binary and real coded genetic
algorithms with six different penalty functions reported in the literature. We also presented a
new method in this paper. Several test problems were solved to demonstrate the use of GAs
in solving BSP. The results were compared and analysed for all seven penalty functions
with binary and real coded GAs. The real coded genetic algorithms worked well compared
to binary coding. Among the penalty functions, the Death Penalty method and Superiority of
Feasible Points were found to be worse. The new method showed its superiority for tight
constrained problems.
9. Acknowledgement
This work is supported by UC special research grant, ADFA, University of New South
Wales, Australia, awarded to Dr. Ruhul Sarker. We like to thank Prof. Xin Yao for many useful
comments during conduct of this research.
10. References
1.
2.
3.
4.
5.
6.
7.
8.
9.
10.
11.
12.
13.
14.
15.
16.
17.
T. Back, Evolutionary Algorithms in Theory and Practice, Oxford University Press, New York, 1996.
A. Cauchy, Method General pour la Resolution des Systeme's d'Equations Simulanees, Compt. Rend., 25
(1847) 536-538.
R. Courant and D. Hilbert, Methods of Mathematical Physic, Interscience, New York, 1 1953.
J. B. Rosen, The Gradient Projection Method for Nonlinear Programming, Part-1: Linear Constraints, J.
Soc. Ind. Appl. Math., 9 (1960) 181-217.
R. Courant, Variational Methods for the Solution of Problems of Equilibrium and Variations, Bull. Am.
Math. Soc., 49 (1960) 1-23.
C. W. Carroll, The Created Response Surface Technique for Optimizing Nonlinear Restrained Systems,
Operations Research, 9 (1961) 169-184.
J. Abadie and J. Carpenter, Generation of the Wolfe Reduced Gradient Method to the Case of Nonlinear
Constraints, In Optimization, R. Fletcher(ed), 1969.
P. Wolfe, Methods of Nonlinear Programming, Notices Am. Math. Soc., 9 (1962) 308.
R. A. Sarker and L. R. Khan, An Optimal Batch Size for a Manufacturing System Operating Under a
Periodic Delivery Policy, Conference of the Asia-Pacific Operational Research Societies (APORS),
Melbourne, 1997.
R. A. Sarker and C. Newton, Determination of Optimal Batch Size for a Manufacturing System, Progress
in Optimization II, Yang, Mees, Fisher & Jennings (editted), Kluwer Academic Publishers (in press),
1999.
A. E. Smith and D. W. Coit, Constraint-Handling Techniques: Penalty Functions, Handbook of
Evolutionary Computation, Release 97/1, IOP Publishing Ltd and Oxford University Press, (1997)
C5.2:1-2:6.
H-P Schwefel, Evolution and Optimum Seeking, Wiley, New York, (1995) 16.
M. S. Bazaraa and C. M. Shetty, Nonlinear Programming: Theory and Algorithms, John Wiley & Sons,
New York, Chapter-9, 1979.
Z. Michalewicz, Genetic Algorithms Numerical Optimization and Constraints, Proc. 6th Int. Conf. on
Genetic Algorithms (Pittsburgh, PA, July 1995) ed L J Eshelman (San Mateo, CA: Morgan Kaufmann),
151-158.
M. Khouja, Z. Michalewicz, and M. Wilmot, The Use of Genetic Algorithms to Solve the Economic Lot
Size Scheduling Problem, European Journal of Operational Research, 110 (1998) 509-524.
R. A. Sarker and C. Newton, Genetic Algorithms for Solving Economic Lot Size Problems, Accepted for
publication in the Proc. of the 26th International conference on Computers and Industrial Engineering,
Melbourne, Australia, December 1999.
A. Homaifar, C. X. Qi, and S. H. Lai, Constrained Optimization Vis Genetic Algorithms, Simulation,
April (1994) 242-253.
13
18. J. A. Joines and C. R. Houck, On the Use of Non-Stationary Penalty Functions to Solve Nonlinear
Constrained Optimization Problems With GAs, Proc. of the IEEE ICEC, (1994) 579-584.
19. Z. Michalewicz and N. Attia, Evolutionary Optimization of Constrained Problems, In A. V. Sebald & L. J.
Fogel (Eds.), Proc. of the 3rd Annual Conf. on Evolutionary Programming, River Edge, NJ: World
Scientific, (1994) 98-108.
20. Z. Michalewicz, Genetic Algorithms + Data Structures = Evolution Programs, 3rd ed., New York,
Springer-Verlag, 1996.
21. Z. Michalewicz, T. Logan, and S. Swaminathan, Evolutionary Operators for Continuous Convex
Parameter Spaces, In A. V. Sebald & L. J. Fogel (Eds.), Proc. of the 3rd Annual Conf. on Evolutionary
Programming, River Edge, NJ: World Scientific, (1994) 84-97.
22. Z. Michalewicz, and M. Schoenauer, Evolutionary Algorithms for Constrained Parameter Optimization
Problems, Evolutionary Computation, 4 (1996) 1-32.
23. D. Powell and M. M. Skolnick, Using Genetic Algorithms in Engineering Design Optimization with
Nonlinear Constraints, In S. Forrest (Ed.), Proc. of the 5th Int. Conf. on Genetic Algorithms, (San Mateo,
CA: Morgan Kaufmann), (1993) 424-430.
14
View publication stats