Genetic algorithms, function optimization, and facility layout design

322

European Journal of Operational Research 63 (1992) 322-346 North-Holland

Genetic algorithms, function optimization, and facility layout design Kar Yan Tam Department of Business Information Systems, School of Business and Management, The Hong Kong University of Science and Technology, Clear Water Bay, Kowloon, Hong Kong Received November 1990; revised December 1991

Abstract: The use of genetic algorithms to solve function optimization problems requires a coding scheme that represents solutions as strings of finite length. In this paper, a coding scheme designed for facility layout is presented. A solution is represented by the postorder sequence of the nodes in a slicing tree. New solutions are generated by applying the various genetic operators in each generation. The genetic algorithm approach is compared with a local search technique in solving problems with sizes ranging from 12 to 30 facilities. The results demonstrate that a genetic algorithm can be a viable tool to solve large scale layout problems. Keywords: Genetic algorithms; facility layout; optimization; Artificial Intelligence

1. Introduction

The problem of facility layout is concerned with the physical placement of a number of interacting facilities on a planar site. It usually arises as a subproblem in the design of an operating system composed of autonomous yet interacting work units. For years, the problem has received considerable attention from firms engaged in manufacturing activities. In setting up a manufacturing system, a facility designer is the person responsible for designing a floor plan layout for the given units. To accomplish the task satisfactorily, the designer has to consider a number of factors, including traffic volume between facilities, area requirements, and geometric constraints of individual facilities. Very often, tradeoffs have to be made between conflicting factors to produce a feasible layout. To support the designer in his task, a suite of computerized layout procedures have been developed (Moore, 1974). A majority of these procedures adopt a problem formulation known as the Quadratic Assignment Problem (QAP). In QAP, the location of each site is specified in advance. The final solution involves assigning one facility to each site so that the total weighted distance between facilities is minimized. Weights can be measured either by an adjacency index or by the volume of materials flow between facilities: the former is usually depicted as a Relationship chart, the latter as a Flow-to chart. It was shown that QAP belongs to the class of NP-complete problems (Sahni and Gonzalez, 1976), and an efficient procedure for solving the problem optimally has not yet been found. During the past three decades, numerous heuristics have been developed to obtain good rather than optimal solutions. In general, these heuristic procedures can be classified into two main categories: construction methods and improvement methods. In construction methods, the output is a sequence representing the order by which facilities are entered into the floor plan. Once entered, a facility is fixed at the allocated location. Correspondence to: K.Y. Tam, Department of Business Information Systems, School of Business and Management, The Hong Kong

University of Science and Technology, Clear Water Bay, Kowloon, Hong Kong. 0377-2217/92/$05.00 © 1992 - Elsevier Science Publishers B.V. All rights reserved

K.Y. Tam / Genetic algorithms, optimization and facility layout design

323

Well known examples of construction procedures are C O R E L A P (Lee and Moore, 1967) and A L D E P (Seehof and Evans, 1967). On the other hand, improvement methods start with an initial layout and attempt to improve it by interchanging facilities. An example of improvement methods is the steepest descent pairwise interchange procedure used by C R A F T (Armour and Buffa, 1963). Compared with construction procedures, systems that use improvement methods tend to produce better quality layouts. Area requirements are not included in the original QAP formulation. To overcome this limitation, systems such as C R A F T include area requirements by dividing the planar site into a grid of square blocks (e.g., 1 m 2) (Tompkins and White, 1984). The space allocated to a facility is represented by a collection of adjacent blocks whose total area equal to that of the facility. To reduce the added complexity introduced by the grid, the number of allowable facility interchanges is largely reduced in these systems. In CRAFT, only facilities having the same area or adjacent facilities sharing a common boundary can be interchanged. This restriction will inevitably limit the search space and will have an adverse impact on the solution quality of the final layout. In addition to the above limitation, existing layout procedures have been criticized for two major problems: (1) The search is conducted from a single point in the solution space, and (2) Geometric constraints are difficult to enforce. While the first problem is related to search techniques in general, the second is more specific to the domain of facility layout. For improvement heuristics, the search is conducted locally from a single point in the solution space. Unfortunately, the search is easily trapped in local optima. A common practice is to execute a number of runs with different initial solutions, then select the best solution as the final result. Simulated annealing (SA) has been suggested to overcome this problem by using a long cooling schedule to approximate asymptotic convergence (Kirkpatrick et al., 1983; Laarhoven and Aarts, 1987; Otten and Ginneken, 1989). The method is motivated by an analogy to a phenomenon in crystallization (Metropolis et al., 1953). Unlike hillclimbing techniques, SA allows occasional down-hill moves to escape from unattractive local optima. It has been applied to solve the location (Skiscim and Golden, 1986) and the traveling salesman problem (Kirkpatrick et al., 1983). The latter has shown to be a special case of the quadratic assignment problem. When compared with a local search method in some facility layout tasks, SA generated better solutions at the expense of prolonged computation time (Tam, 1991). In each iteration, both techniques search only a single region defined by the neighborhood of an existing solution. In SA, the attempt to sample different regions of the solution space may only be realized with a cooling schedule that has (1) a high initial temperature, (2) a small temperature decrement in each iteration, and (3) a long annealing chain. All these combine to form a schedule that will demand a prohibitively long CPU time. Most importantly, both SA and local search methods do not learn from the solution space already explored; consequently, useful information about the function surface that can be inferred from known solutions are left unused. The second problem is related to the model formulation itself. The inability to describe and enforce geometric constraints has been criticized as a major deficiency of existing systems. Overlaying a grid on a floor plan only provides a means to satisfy the area requirements, it does not specify or enforce any shape constraints commonly found in a layout task. Shapes that are infeasible for some facilities may be generated, resulting in a layout of little practical use. In this paper, we address these two problems by using a symbolic representation of a floor plan layout and outlining a solution procedure using genetic algorithms. The various factors mentioned earlier are all taken into consideration in our formulation. The objective is to minimize the weighted flow of traffic while satisfying the various area and geometric constraints of individual facilities. The use of genetic algorithms allows an effective search of the solution space by sampling different regions of the space in parallel. In Section 2, an overview of genetic algorithms is presented. This is followed by a discussion on the layout representation used in our problem formulation in Section 3. In Section 4, computational results of the procedure are reported. Section 5 concludes the paper with a discussion on this method's limitations, extensions, and potential for future work.

K.Y. Tam / Genetic algorithms, optimization and facility layout design

324

2. Genetic algorithms Genetic Algorithms (GA) is a class of iterative procedures that simulate the evolution process of a population of structures subject to the competitive forces prescribed in Darwin's 'survival of the fittest' principle. The process of evolution is random yet guided by a selection mechanism based on the fitness of individual structures. If carefully designed and properly implemented, a G A will exhibit a behavior similar to that described in Darwin's evolution theory - relatively high fitness structures have a larger chance to survive and to produce even higher fitness offsprings. The result will be an increase in the overall fitness of a population in each new generation. It seems appropriate to draw an analogy between the evolution process simulated by a G A and the search process conducted in solving a function optimization problem. Both involve searching through a set of alternatives. The link between a G A and a search procedure is shown below: structure o fitness o

solution, objective function value,

evolution ~, search. Suppose a one-to-one mapping between the structure and the solution space can be found; a point in the solution space can be coded as a structure in a GA. Given a proper transformation, the fitness of a structure can be interpreted as a proxy of the original objective function. While most search techniques follow a point search mechanism, a G A maintains a population of points to be explored. In traditional search methods such as hillclimbing, the search process is originated from a single point (i.e., a neighborhood search) in each iteration. Little or nothing is learned during the search process even though useful information about the landscape of the function can be uncovered from the already explored search space. A G A goes beyond a mere sequential search of each point (structure) in the population by identifying and exploiting common building blocks of good structures in the population. These building blocks correspond to regions in the space where good solutions are likely to be found. A G A uses a set of genetic operators to select, recombine and alter existing structures according to sound principles to direct the search towards these regions. Holland (1973) has shown the striking similarity between G A and the well known k-bandit problem in the way their search processes are directed. Both problems involve deciding the amount of search efforts allocated to different regions of the solution space as the search proceeds. While G A has been used in numerous areas (Brooker et al., 1989; DeJong, 1988), we will focus on its application in function optimization. A general genetic algorithm is shown in Figure 1.

GA() {A genetic algorithm to minimize F} begin t =0; initialize Pt; evaluate Pt; While not (terminate condition) do begin t=t+l; select Pt from Pt-1; {reproduction operator} {crossover and mutation operator} recombine Pt; evaluate P. end end. Figure 1. A genetic algorithm

K. E Tam / Genetic algorithms, optimization and facility layout design

325

Implementation of a G A presupposes the following: (1) A mapping between the solution space and the structure space; (2) A fitness measure of structures; (3) A set of genetic operators. In a GA, a structure s is coded as a string of length k (k > 0) over an alphabet set V. The structure space is defined as the set of strings belonging to V k. Thus, the size of the structure space is I V [ k where ] V I is the n u m b e r of symbols in V. A popular choice for V is the binary set {0, 1}. A crucial element in the study of G A is the concept of schema. A schema is a string of length k defined over the alphabet set V U {#}. The symbol # is a 'don't care' or wild card symbol that can be used as a substitute for any symbol in V. Since a schema corresponds to a plane in the hypercube defined by the Cartesian product of V, it is also called hyperplane in G A literature. For instance, given a structure s = 011011 and a schema H = 0#1011, we say that s is a representative of H because s can be derived from H by substituting the # symbol in H for '1'. It follows that the size of the schema space is (I V I + 1) k. The number of representatives of H is denoted by M(H). In our previous example, M ( 0 # 1 0 1 I ) is two. Two other important properties of a schema H are defined below: (1) O r d e r of H, denoted by o ( H ) : the number of fixed symbols in H (e.g., o(0#1011) = 5). (2) Defining Length of H, denoted by 6 ( H ) : the difference between the first and the last fixed symbol in H (e.g., 3 ( # 1 # 0 1 # # ) = 5 - 2 = 3). The fitness of a structure is measured by a function p. defined as /z : S -~ R + where S is the set of all structures (i.e., V k) and R + is the set of non-negative real numbers. If it is known that the function f of the underlying optimization problem is always positive, then f can be directly used as/z. Otherwise,/z is a transformation of f. The actual form of transformation depends on a number of factors: (1) whether max f or min f is used, (2) the selection mechanism employed, and (3) the scaling function used. Detailed discussions of these factors are covered in Goldberg (1989) and will not be repeated here. When applied to solve function optimization problems, a G A works as follows: A population of solutions coded as strings of fixed length are maintained during the search. In each iteration, a new population P,+t is created by retaining old solutions and generating new solutions from the previous population Pt. This is accomplished by applying genetic operators to the solutions in /'i. Three basic operators form the core of most G A implementations: (1) reproduction operator, (2) crossover operator, and (3) mutation operator. O t h e r operators have also been suggested but they are either derivatives of these core operators or specific operators targeted for a particular problem. The overall quality of solutions will improve after each iteration. The iterating process stops when a certain level of solution quality is reached or convergence is observed. Reproduction operator: Under the reproduction operator, the chance of being selected to remain in the new population P,+t is proportional to a structure's fitness. This operator assigns to each structure a sampling rate T(s, t), defined as the expected number of offsprings to be generated from that structure in generation t. To coincide with evolution theory, given two structures s' and s", if iz(s') > #(s"), then T(s', t) > T(s", t). The most widely adopted definition of the sampling rate is T(s, t) = ~ ( s ) / ~ ( s ) where the numerator is the fitness of s and the denominator is the average fitness of population P,. Notice that structures with above average fitness will have a higher survival probability than those with below average fitness. Crossover operator: With only the reproduction operator, the population will become more and more homogeneous after each generation. The crossover operator is included in a G A for two purposes. First, it introduces new structures by recombining existing ones. Second, it serves as a sieve to eliminate low fitness schemas. It operates as follow: given two structures s' and s", they interchange a substring of themselves according to a crossover point to form two new structures. Suppose s' is 011 I 11 and s" is 101 r00 with a crossover point location indicated by I. By interchanging the substring to the right of the crossover point, two new structures 01100 and 10111 are created. To avoid chaotic behavior, not all structures in the new population are generated by the crossover operator. The probability of applying this operator, or simply the crossover rate, is denoted by PcMutation operator: This operator introduces random changes to the structures in a population by

K.Y Tam / Geneticalgorithms, optimizationandfacility layoutdesign

326

changing a symbol in a structure with a probability (or mutation rate) Pm" For example, if Pm is 0.01, then in each generation, there is a 1% change that a structure in the population will be changed by alternating one of its symbol. The effect of this operator becomes more apparent when the crossover rate is small. In this case, the mutation operator reduces the possibility of premature convergence. In addition to the above items, a complete G A implementation involves specifying a number of parameters, including population size, number of generations, scaling function, generation gap and selection strategy. The power of G A can be explained by two major theoretical results. First, Holland's (1975) schema theorem states that the expected number of representatives of a schema H in Pt+l, denoted by M(H, t + 1), using the three operators described above is given by the following expression:

M(H,t+I)=M(H,t)

I~(H) 1 - P c

~-1

o(H)p m •

The theorem states that H will receive an increasing number of representatives in Pt +1 if it has high fitness (/~(H)), short length (6(H)), and low order (o(H)). Since the exact value o f / z ( H ) is unknown, it is estimated by the average fitness of H ' s representatives in ft. Grefenstette and Baker (1989) challenge the appropriateness of the schema theorem in dealing with fitness functions that are nonlinear transformations of the original function. They argue that the reproduction operator described above is sensitive to the transformation used. Nevertheless, Holland's results still hold for most implementations that use a fitness function of the f o r m / ~ ( x ) = af(x) + b where a and b are constants with b > 0 if max f and b < 0 i f m i n f. The second result is what Holland calls the intrinsic parallelism property of GA (Holland, 1980). According to the schema theorem, each schema with at least one representative in Pt will receive M(H, t + 1) representatives in Pt+l. The result states that the number of schemas in a population is at least O(n 3) where n is the number of structures in a population. In other words, a G A allocates search effort to O(n 3) distinct hyperplanes (i.e. schemas) in the solution space in parallel. Most importantly, this search process is conducted implicitly. The ability to explore different regions of a solution space in parallel is a unique feature of GA, and the solution quality obtained may justify the overhead of maintaining a population of solutions.

3. Problem formulation

In this section, we will describe the details of the proposed problem formulation.

3.1. Symbolic coding of a floor plan layout A challenging problem in applying G A is the coding of solutions as strings of finite length. Care must be exercised to ensure that the meaning of the original problem is preserved under the coding scheme. Goldberg (1989) arrives at two principles for designing G A coding: Principle of meaningful building blocks: " T h e users should select a coding so that short low-order schemata are relevant to the underlying problem and relatively unrelated to schemata over the other fixed positions". Principle of minimal alphabets: " T h e user should select the smallest alphabet set that permits a natural expression of the problem". A floor plan layout is represented as a slicing structure constructed by recursively partitioning a rectangular block (i.e., the floor plan) in such a way that each rectangular partition in the slicing structure corresponds to the space allocated to a facility. Table 1 shows the geometric constraints of 6 facilities. Column 2 in the table contains the area requirement of each facility. Columns 3-5 will be discussed later. A layout of these 6 facilities is represented as a slicing structure shown in Figure 2. We

K.Y. Tam / Genetic algorithms, optimization and facility layout design

327

Table 1 Geometric constraints of 6 facilities Facility id

Area

1 2 3 4 5 6

Deadspace ratio

Aspect ratio

100 80 50 60 120 40

lower bound

upper bound

upper bound

0.5 0.4 0.3 0.1 0.9 0.4

1.2 1.1 3 5.5 1.9 1.2

0 0 0 0 0 0

say the area requirement is satisfied if the area allocated to each facility is at least as large as the area required. Notice that the area requirements are all satisfied in Figure 2. An equivalent representation of a slicing structure is a slicing tree. A slicing tree is a binary tree which shows the recursive partitioning process that generates a slicing structure. An example of a slicing tree is shown in Figure 3. Each internal node represents the way a rectangular partition is cut. The kind of cut is denoted by a letter assigned to each internal node (e.g., 'u', 'b'). As shown in Figure 3, partitions reserved for facilities reside in the leaves of the tree. Each leaf is assigned an unique integer corresponding to the identifier (id) of a facility. Suppose we fix the structure of a tree and only change the cut of an internal node, a different layout is created. For example, by changing the internal nodes one at a time, we obtain a sequence of layouts as shown in Figure 4. The space of all layouts is defined as the set of all slicing trees S that can be generated by rearranging the cuts of a slicing tree. Notice that each element s e S has the same tree structure. They only differ in the cuts of their internal nodes. Given a tree structure, a more convenient representation of a slicing tree is the postorder sequence of its nodes. For example, the slicing tree shown in Figure 3 can be represented by the string 65u4132ulbr. It is quite clear that two trees are identical if and only if their structures and their postorder sequences are the same. Since this representation is similar to the way operators and operands are stated in postfix arithmetic expressions (e.g., 1 2 . 3 + ), we will refer to the cut symbols as operators and the facility ids as operands. Four cut operators are used in the procedure: l for left cut, r for right cut, u for up cut, and b for bottom cut. Each operator is explained pictorially in Figure 5. Although it seems redundant to have two operators instead of one for each dimension, this is necessary in order to maintain the order of the operands in a

12.78

8.89 6

7.82

1

3.91

3

3.33 I..5

18 5

6.25

13.5

2

25 Figure 2. A slicing structure

Figure 3. A slicing tree

328

K.Y. Tam / Genetic algorithms, optimization and facility layout design

4

%

3

2

6 1

5

2

4

Figure 4. A sequence of slicing trees

sequence. For instance, consider changing the internal node with label u to b in Figure 6. If only one operator is used for each dimension, say v for the vertical axis and h for the horizontal axis, the node's children have to be interchanged in order to reflect this change. This in turn will alter the order of operands in the sequence. Using two opposite operators along each dimension (e.g., l and r, u and b), the operands are fixed in the sequence. Since the operands remain unchanged, we can delete them and retain only the operators. For example, the sequence 6 5 u 4 1 3 2 u l b r can simply be stated as ulubr. Given a tree structure, the following equivalent representations are observed: slicing tree

¢* operator sequence

and floor plan ~

slicing structure.

~

structure in G A

K.Y.. Tam / Genetic algorithms, optimization and facility layout design

a

2

329

1

1

2

£L

2 1

1

2 Figure 5.

Figure 6. Cut tree with 2 operators

Using the terminology of GA, a structure is the same as an operator sequence. All structures are generated from the alphabet set {b, u, r, l} n - I where n is the number of facilities. The size of the solution space is therefore 4 n-1. Each alphabet (i.e. operator) in a structure represents a unique way to cut a partition and is relatively unrelated to alphabets in other positions (principle of meaningful building blocks). It is also shown above that the set is minimal without distorting the structure of a slicing tree (principle of minimal alphabets). 3.2. Constructing a layout from an operator sequence To complete our discussion of layout coding, a procedure to convert a slicing tree to a slicing structure is outlined. Since the procedure must allocate enough space to each facility, we will assume that the total usable floor area is at least as large as the total area of the facilities combined. To enforce the area requirements, the location where a rectangular partition is cut, called the cut point, must be decided so that the split partitions receive the area desired. The cut point is determined in a t o p - d o w n fashion starting from the root of the tree and passing down to the children nodes recursively. Since the dimensions of a node (i.e., a partition) and the areas of its children are known, the cut point can be located in a straightforward manner, provided the partition does not contain any occupied areas. The latter are regions on the floor plan that have already been allocated for other purposes. The cut point

K.Y. Tam / Genetic algorithms, optimization and facility layout design

330

must be searched for in situations where a partition overlaps with some occupied areas. Binary search is used to locate the cut point in these cases. The search can be performed at any arbitrary level of precision as specified by the facility designer. An upper bound on the time complexity of finding all cut points of a slicing tree (operator sequence) with n internal nodes (operators) is O(nm log m) where m = max{wf, hf}/y, where wf and hf are the width and height of the floor plan, respectively, and 3' is a parameter governing the precision of the search. If the floor plan is clear of occupied areas, then the time complexity is reduced to O(n).

3.3. Slicing tree structure construction Our formulation is based on a fixed slicing tree structure. The structure must be constructed in a way that reflects the proximity requirements of interacting facilities. The approach used here involves applying numerical clustering techniques to place facilities with large traffic volumes closer together. The use of clustering techniques requires a distance (or dissimilarity) measure between facilities. A distance measure, denoted by Agj, between every pair of facilities is constructed using the traffic information 1

Aij= l + (vij + vji) where u~ and vii are the traffic flows from facility i to j and from j to i, respectively. For every pair of facilities i and j, their distance Aij is stored in a symmetric distance matrix. This matrix is input to a numerical clustering procedure to generate a dendogram. The dendogram depicts the hierarchical process of clusters merging and can be used to provide the structure of a slicing tree. Various clustering procedures, such as single linkage, complete linkage, and density linkage, have been developed (Anderberg, 1973) and are easily assessable in a number of statistical packages (e.g., SAS, SPSS). Although the single linkage method has been used extensively in cellular systems design (Chu, 1989), we found it less attractive than the average linkage method in constructing a slicing tree structure. Single linkage method tends to generate trees that are deeply skewed, especially in cases where some elements in the distance matrix are identical. Using this method, two clusters are merged if they have the closest pair of facilities. Ties are resolved arbitrarily. It may happen that a cluster is merged with a sequence of singleton clusters if these clusters are the same distance from some elements of the initial cluster. To illustrate this, consider the traffic volume matrix of fifteen facilities shown in Table 2. If we convert the dendogram constructed using the single linkage method into a slicing tree structure, a very skewed tree is formed (Figure 7). A skewed tree should be avoided because it is quite likely that very thin partitions will be formed thus violating the geometric constraints of these facilities. Using the average linkage method, one obtains a more balanced tree structure (Figure 8). In this method, the distance between two clusters is the average distance between pairs of facilities, one in each cluster. Average linkage tends to join clusters with small variance and avoids the skewness introduced by the single linkage method. Another clustering method, called Ward's minimum variance method (Ward, 1963), gives similar results. The use of a dendogram to provide a structure for the slicing tree allows facilities with large traffic volumes to be placed closer to each other. Note that facilities having the same parent in the slicing tree will be placed adjacent to each other, irrespective of the kind of cut used.

3. 4. Geometric properties of facilities The shape of a facility is described by two parameters: aspect ratio and deadspace ratio. The aspect ratio a i of facility i is defined as height of the partition allocated to facility i a; = width of the partition allocated to facility i

K.Y. Tam / Genetic algorithms, optim&ation and facility layout design

331

Table 2 S y m m e t r i c t r a f f i c m a t r i x o f 15 f a c i l i t i e s F a c i l i t y id

Facility id 1

2

3

5

6

7

8

9

10

13

14

15

5

1

-

1

2

2

2

2

-

4

-

-

1

3

2

2

2

3

2

-

2

-

10

5

-

-

10

2

-

2

5

4

5

2

2

5

5

5

-

1

1

5

-

-

2

1

-

2

5

-

-

3

5

5

5

1

-

3

-

5

5

-

2

2

1

5

-

-

2

5

10

-

6

-

1

5

5

5

1

-

-

5

2

10

-

5

-

-

9

-

10

5

10

-

2

10

-

-

4

-

-

5

-

5

-

5

-

-

3

3

-

10

2

-

4

1

10

2

-

3 4 5 6 7

4

8

11

11

12

12

13 14 15

Figure

8. S l i c i n g t r e e

constructed

by t h e

average

linkage

method

F i g u r e 7. S l i c i n g t r e e c o n s t r u c t e d by t h e s i n g l e l i n k a g e m e t h o d

O n e can impose constraints on the shape of a facility by restricting its aspect ratio t o [ a mi i n , a ml a x ] , where i n and area i x are the lower and u p p e r b o u n d of the aspect ratio a i. W e can also impose constraints on ami facility orientation by classifying facilities into two categories: free oriented facilities and fixed oriented facilities. A free oriented facility allows both vertical and horizontal orientation. For these facilities, their aspect ratios can either be a i or 1 / / a i. Thus, the feasible aspect ratio range b e c o m e s

aspect ratio range of free oriented facility i =

min

aiin,

i amax

,

max

a ....

i

amin

•

K.Y. Tam / Genetic algorithms, optimization and facility layout design

332

::...'.:.3!~•

"

2

10 3

I

1 3

7 7

4

3

F i g u r e 9. F l o o r p l a n w i t h o c c u p i e d s p a c e s

A fixed oriented facility allows only vertical (or horizontal) orientation, and its allowable aspect ratio range is simply i i aspect ratio of fixed oriented facility i = [amin, amax].

The range of aspect ratios represent the constraints imposed by the type of work performed in each facility. In a manufacturing environment, each facility may represent a manufacturing cell with its shape constraints defined by the way machines are arranged in the cell (Heragu and Kusiak, 1988). Consider the floor plan in Figure 2. The aspect ratios are 0.61, 0.49, 0.31, 5.41, 1.52, and 0.51 for facilities 1-6, respectively. Note that all aspect ratios are within the range specified in Table 1. Very often, the entire floor plan cannot be used for facilities placement. There are always locations on the floor plan that have already been occupied by existing facilities such as elevator shafts, utilities and columns, etc. These occupied spaces must be subtracted from the area allocated to a facility. Occupied spaces are modelled as one or more adjacent rectangular blocks. Consider the three occupied spaces shown in Figure 9: the area allocated to facility 1 (drawn in dark perimeter) equals the area of the rectangular partition minus the portion of the occupied space overlapping it. In general, all overlapping occupied spaces must be subtracted from the partition allocated to a facility. Occupied spaces can also be used to mold the perimeter of the floor plan to approximate any shape. This is accomplished by placing 'arbitrary' occupied spaces on the boundary of a rectangular block. The existence of occupied areas inside a facility will distort the rectangular shape of the usable area, therefore it is imperative to provide a measure which reflects the degree of shape distortion. The deadspace ratio is introduced for this purpose. The deadspace ratio of a facility i denoted by o i is defined as area of total occupied spaces in the partition allocated to facility i 0i

area of the partition allocated to facility i

If the floor plan is free of any occupied areas, the deadspace ratios of all facilities will be zero. For each facility i, one can restrict the value of o i to fall within [0, Oma x i ] where Oma x i is the upper bound of o i. Consider the boundary of facility 1 in Figure 9. The area of the rectangular partition allocated to 1 is 10- 7 = 70. The portion of the occupied space that overlaps with facility 1 has an area of 9, therefore its deadspace ratio is 9 / 7 0 ~ 0.13. Based on the above formulation, given a slicing tree, one can compute (1) the boundaries, (2) the aspect ratios, and (3) the deadspace ratios of all facilities. In fact, all these can be computed once the cut points of a tree are known.

K. Y. Tam / Genetic algorithms, optimization and facility layout design

333

3.5. P r o b l e m d e f i n i t i o n

The objective is to allocate space to facilities in such a way that facilities with large traffic volumes are placed on close proximity while still satisfying the area and shape constraints of individual facilities. Given a slicing tree s ~ S, the space allocation problem is formulated as

minF= s~S

~ ~vijdij i--I j ~ l

subject to

i i a m i n ~< a i <_ a m a x , i 0 ~ 0 i <_~Omax,

i,j=l,2

. . . . ,n

where = = = = = = =

Traffic volume between facilities i and j. Rectilinear distance between the centroids of the partitions allocated to facilities i and j. dij Aspect ratio of the partition allocated to facility i. ai Deadspace ratio of the partition allocated to facility i. oi i Lower bound of a i. amin i U p p e r bound of a i. a max l U p p e r bound of o r Omax G A s are designed to solve unconstrained optimization problems. In our formulation, the above constraints are converted into a penalty function. The problem is reformulated as an unconstrained optimization problem stated as vij

min s~S

F= ~ i-I

where

OLk = m a x J[~k =

~ vijdij q- ~

j=l

(w~olk -~-w~k~k)

k=l

( (( O, m a x

a k

-

-

max akax,

1

-

k a min -

})(

, min a in,

1}

k

--ak

)))

,

amax

k

max{0, O k -- Ornax},

wL w~ >_0, i,j,k=l,2

..... n.

The first term of F measures the weighted traffic volume. It is defined as the product of the rectilinear distance dij between two facilities and their traffic volume vii. The second term represents a penalty function for geometric constraints. T h e first term of the penalty function measures the extent to which the aspect ratio constraints are violated. Since each facility has a unique tolerance level of violating this constraint, a k is multiplied by a weight factor w~ >_ 0. For fixed oriented facilities, the aspect ratio penalty ak can be further simplified to a k = max{O, max{a k - amax, k a mk i n - - a k } } • Similarly, the second term flk penalizes the existence of occupied spaces within a facility. Individual weights w~ > 0 are assigned to facilities to reflect their tolerance levels for irregular shapes.

K.Y Tam / Genetic algorithms, optimization and faciflty layout design

334

4. Implementation

The genetic algorithm was coded in Pascal and implemented on an EMX machine. Four layouts with facilities ranging in number from 12 to 30 were retrieved from Nugent et al. (1968). Since traffic volumes are symmetrical, only distinct pairs of facilities are considered (i.e., ui~dij, i = 1 , . . . , n, i
•

(area of facility i ) = total usable area of the floor plan.

i=1

In cases where there are more usable areas than required, a space amortization technique is used. That is, the extra space is amortized to each facility in proportion to its area. As shown in the figures, there are two occupied spaces in the 12-facility, three in the 15- and 30-facility, and four in the 20-facility layout. The slicing tree structure of each layout was constructed by the average linkage method and are shown in Figures 14-17.

Table 3 Area and Aspect ratio of each facility Facility id

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

Area

100 80 50 60 120 40 20 40 150 120 50 10 20 30 50 20 40 20 80 100 40 50 80 10 40 10 40 10 80 40

Aspect ratio lower bound

upper bound

0.7 1 0.7 0.5 0.9 0.6 0.7 1 0.8 0.5 0.7 0.8 0.95 0.75 0.9 0.8 0.4 0.9 1 0.95 0.5 1 0.6 0.9 0.8 0.5 0.8 0.5 0.9 0.9

1 1 1.3 0.8 1 1 1.4 1 1.1 1.5 1.1 1.2 1.5 1.25 1.1 1.5 1.4 1.9 1 1.15 1.5 1.1 1 1 1. I 1.2 1 1.3 1.05 1.1

K. Z Tam / Genetic algorithms, optimization and facility layout design

335

35

81

10

8

5

17

32

Figure 10. Floor plan for 12 facilities

Values of the parameters used in the experiments are: population size N = 30, crossover rate Pc = 0.95, and mutation rate Pm = 0.01. The selection of parameters is based on two extensive empirical studies by Grefenstette (1986) and Schaffer et al. (1989). The search for an optimal p a r a m e t e r set has been a major issue in the study of genetic algorithms. Initial work in this area was taken by DeJong in the mid 70s. In his seminal dissertation, DeJong (1975) compared different combinations of N, Pc, Pm, and a number of other parameters in a test bed of functions. Two performance measures were introduced to gauge the effectiveness of a search strategy. According to DeJong, the on-line performance of a search strategy is defined as the average fitness of a population up to T where T is the number of trials attempted to that point, while off-line performance is defined as the average best fitness of a population up to T. Empirical results generated from his study indicate that the following p a r a m e t e r set yields the best on-line performance: N = 50-60, Pc = 0.6, and P m = 0.001. In an effort to find the optimal p a r a m e t e r set in solving a metalevel G A problem, Grefenstette (1986) performed a number of experiments which include other parameters such as generation gap, scaling and selection strategy for a m o r e extensive comparison. The following values were reported to generate the best on-line performance among all combinations studied: N = 30, Pc = 0.95, and Pm = 0.01.

:54

10

5 4

I0

25

15

lO

5l

°I 10

34

Figure 11. Floor plan for 15 facilities

/~Y. Tam / Genetic algorithms, optimization and facility layout design

336

20

10

10

10

10 25 30 20

I

10

$0 Figure 12. Floor plan for 20 facilities

Grefenstette's results were reinforced by Schaffer et al. (1989) in a factorial experiment designed to study the effects of different parameter combinations on on-line performance. Four parameters with a total of 840 cells were included in the study. The results suggest the following values offer the best on-line performance: N = 20-30, Pc = 0.95, Pm = 0.01.

35

5 I0

30 25

10

105

10

30

Figure 13. Floor plan for 30 facilities

Figure 14. Slicing tree for 12 facilities

K.Y. Tam / Genetic algorithms, optimization and facility layout design

Figure 15. Slicing tree for 15 facilities

337

Figure 16. Slicing tree for 20 facilities

Since there is supporting evidence that the reported p a r a m e t e r sets are function independent, they are also applicable to the objective function described in Section 4. In addition to the suggested values for N, Pc and Pm, we adopted the common practice of using scaling function (linear scaling is used with a scaling factor of 1.8) to avoid p r e m a t u r e convergence. A generation gap of 1 was also used. The results of G A were c o m p a r e d with that of a hillclimbing method denoted by HC. H C searches through a neighborhood defined as the set of operator sequences derived from changing an operator. To

Figure 17. Slicing tree for 30 facilities

338

K.Y. Tam / Genetic algorithms, optimization and facility layout design

Table 4 Results of 12-facilitylayout No.

1 2 3 4 5 6 7 8 9 10 min max avg std

Genetic algorithm

Hillclimbing

Cost

No. of tran.

Sequence

Cost

No. of tran.

Sequence

6625.63 6926.06 6729.16 6935.75 6719.75 6896.54 6390.34 6719.75 6702.09 6607.68 6390.34 6935.75 6725.28 166.833

4500 4500 4500 4500 4500 4500 4500 4500 4500 4500

lulrlbubblb llurbbruubl lubrlrrbblb lllruuubbbu lbbrlrruulu rublrulbbrl rurlrbubbbb lbbrlrruulu ulrrbllrrrr brllbrrurll

6896.42 6839.55 6957.69 6943.91 6968.97 7167.05 7007.07 6589.15 6948.31 7275.14 6589.15 7275.14 6959.33 182.86

4541 4510 4564 4546 4619 4537 4553 4518 4525 4557

llblrulbbrl urbrlulbbrr blrlbrlrrrr urrrrlrbrbr lbbrlbbubbb ulurubbbbbb bubrlblbulr rbrbblrrrrr bbbubrrlrll llrburrlull

provide a fair c o m p a r i s o n , the n u m b e r of solutions g e n e r a t e d by H C was limited to 4500 or to the closest local o p t i m a l solution, whichever was larger. T o d e m o n s t r a t e the effectiveness of G A as a g e n e r a l search t e c h n i q u e , we did not include specific g e n e t i c o p e r a t o r s tailored for the layout r e p r e s e n t a t i o n used in our study. By simply a d o p t i n g a p a r a m e t e r set suggested in the literature, we have avoided any exploratory e x p e r i m e n t s which would inevitably increase the o v e r h e a d of o u r i m p l e m e n t a t i o n . A r a n d o m initialization process was used in o u r i m p l e m e n t a t i o n in o r d e r to avoid any biases that might affect the comparative study. Initial solutions were o b t a i n e d by r a n d o m l y g e n e r a t i n g cut o p e r a t o r s for 30 slicing trees. Note that o n e could improve the overall quality of the initial p o p u l a t i o n by g e n e r a t i n g solutions in a m o r e intelligent way, O n e alternative is to have the p o p u l a t i o n initialized by a local search t e c h n i q u e .

Table 5 Results of 15-facilities layout No.

1 2 3 4 5 6 7 8 9 10 min max avg std

Genetic algorithm

Hillclimbing

Cost

No. of tran.

Sequence

Cost

No. of tran.

Sequence

10562.67 10401.8 10156.41 10582.82 10378.14 10287.54 10144.55 10308.03 10099.08 10404.39 10099.08 10582.82 10332.54 167.04

4500 4500 4500 4500 4500 4500 4500 4500 4500 4500

rbrrrbbblrlbur ruulllrlllburr rbrruubllrlbbr lulbulrubublru ruubrrlrrrbull rbbrrrlrrrubl rurrbblllrbuur rbrrbrlbuburrb lulluuburlllrl rurruubllrlubr

11156.82 11879.45 11175.29 11470.74 10751.45 11135.97 10908.02 11168.30 11082.09 11200.62 10751.45 11879.45 11192.88 305.94

4588 4538 4555 4524 4501 4553 4518 4529 4539 4553

lbrblbrllllbuu bbruublrllbbul Illbbrbruuluub ubrubrrlbubblu lbbbbbubbrllll lbrubbublrrulr bblubblbrrulU ruullburrrblbl uulrbbrrrrbbuu bbbulbrulrulrr

K.Y. Tam / Genetic algorithms, optimization and facility layout design

339

8000

8 E

.e

7000

No. of generations 6000

!

I

100

20O

Figure 18. Minimum cost vs. generations (12 facilities)

Each G A was run for 150 generations with 10 different sets of initial solutions. Statistics of each run were gathered. These included the best solution and the average solution in each generation. Both performance measures were applied to the observed scaled fitness and the actual objective function. The results of the four layouts are shown in Tables 4-7. G A outperformed HC both in terms of minimum and average costs. As shown in the tables, the observed costs of G A were less random than HC in all layout tasks. The standard deviation of G A in each experiment was lower than that of HC. In the 15- and 30-facility layout, HC's standard deviation

12000

\ 8 11000

No. of generations 10000

i

!

100

2OO

Figure 19. Minimum cost vs. generations (15 facilities)

K.Y. Tam / Genetic algorithms, optim&ation and facility layout design

340 30000

29000

8 E 28000

.E

27000

No. of generations 26000 0

!

I

100

200

Figure 20. Minimum cost vs. generations (20 facilities) 58000

56000

54000

8 E •~

s2ooo

50000

48000

No. of generations 46000

]

I

0

100 Figure 21. Minimum cost vs. generations (30 facilities)

lO 2 11

6

12 I

3

Figure 22. Floor plan layout (12 facilities)

200

342

K.Y. Tam / Genetic algorithms, optimization and facility layout design

10

6I

I1

2

8

6

I7

15

19

13 2O

I 5

4

Figure 24. Floor plan layout (20 facilities)

was twice that of GA's. When the number of facilities is small (n = 12), the improvement is only marginal. However, as n increases, noticeable improvements can be observed. Percentage reduction in minimum costs of 10.5% and average costs of 13% were achieved when n = 30. The best cost vs. time (measured in generations) of the best run in each experiment are shown in Figures 18-21. In the 12- and 15-facility layout, the cost remained virtually the same after 100 generations. Similar converging behavior was observed in the 20-facility layout. When compared with HC, G A utilized only two-thirds of its computation time in the first three experiments. Even so, G A still outperformed HC. In the 30-facility layout, convergence of the best cost is less apparent, and GA scored the largest improvements over HC in this case.

0

I

1261

10

21 29

6

I 13

22

16

7

I 8I

25

28

19

3

,.

23

24112 17 15

11

4

14

20

27 Figure 25. Floor plan layout (30 facilities)

I

K.E Tam / Genetic algorithms, optimization and facility layout design

343

Notice that such an improvement was achieved, even HC attempted more solutions than GA. Furthermore, the parameters set used in our G A implementation has not been fine-tuned to fit the symbolic layout representation. We would expect further improvement by developing tailored genetic operators and tuning the parameter values. The best layouts obtained by GA are shown in Figures 22-25.

5. Discussion and conclusion

We have presented a genetic approach to the facility layout problem in this paper. The symbolic layout representation takes into account both the area and shape constraints of individual facilities. Unlike most existing layout procedures where the floor plan is divided into a grid, facilities are laid on a continual plane. Direct comparisons with existing systems such as CRAFT, C O R E L A P , and A L D E P are difficult because their objective functions are different. According to Tompkins and White (1984), none of these systems allows a facility designer to specify constraints on the shape of individual facilities. C O R E L A P and A L D E P are constructive systems which take in facilities one at a time. The final shape of each facility in A L D E P is determined by the width of the vertical scanning path and the order by which facilities are placed on the site. If unusual shape borders are identified in the final layout, it is discarded and a new layout is attempted by changing either the path width or the input sequence. C O R E L A P uses placing rating and boundary length to determine the location of a facility. The latter measures the length of the boundaries common to the new facility and all facilities already in the layout. One may specify a global constraint on boundary length, but no explicit control can be exercised on the shape of individual facilities. Furthermore, the boundaries are very sensitive to the layout scale and the placement sequence. C R A F T uses an improvement approach by performing pairwise interchanges. Interchanges can be made only when two facilities have the same area or they share the same boundary. Since there is no control over the sequence of interchanges, infeasible layouts may be generated. The final layout is also sensitive to the 'resolution' of the floor plan grid. Like C O R E L A P , while it is possible to specify a global boundary constraint, such a practice will further restrict the search space of CRAFT. The layout representation presented in this paper guarantee rectangular shapes, provided the floor plan is free of occupied spaces. The recursive partitioning mechanism will generate layouts to any degree of precision, thus avoiding the scaling problem. The aspect and deadspace ratios are unique to each facility, and shape violation tolerance can be represented by its associated weights (i.e., wf, and w~). The four operators restrict the solution space to layouts that are generated by recursively dividing the floor plan into rectangular partitions. Obviously, this set of operators is by no means complete. Yet, it illustrates a general representation that can be extended to include additional shapes. For instance, we can define 4 new operators that handle L-shape facilities (Figure 26) as follows: LNE: The enclosing facility is located at the north-east corner. LNw: The enclosing facility is located at the north-west corner. LsE: The enclosing facility is located at the south-east corner. Lsw: The enclosing facility is located at the south-west corner. New operators may introduce additional constraints. Some of them can be incorporated in the penalty function while others must be enforced in the way solutions are generated. Although direct comparisons with existing systems are difficult, it is worthwhile to illustrate the different solution generation mechanism in system like C R A F T and in the proposed procedure. Consider laying out three facilities with area 8, 12, and 10 units, respectively, on a 3-by-10 floor plan. Given an initial solution, C R A F T will generate three solutions in the first iteration as shown in Figure 27. After that, only two solutions will be generated in subsequent iterations. On the other hand, the proposed procedure is able to explore any point in the entire solution space irrespective of the initial population of solutions. Notice in Figure 28 that all facilities are rectangular in shape as compared with

K. Y Tam / Genetic algorithms, optimization and facility layout design

345

:::::::::::::::::::::::::::::::::::::::: ~,.\\\\\\\\\\\\\\\\\\\\\\\\\\\\\~ I 32u lu

32ulb

3211u

32flu

I 3261u

3211b

3~rlh

32blb

3211r

32xlr

32blr

3~11

32bll

~\\\\\\\'~\\\\\\\~\\\~ ,~'-..\\\\\\\\\\\\\\\\~ 32ulr

32ull

32111

Figure 28. Solution space of the 3-facility layout problem

ness for optimizing off-line performance may be quite different. Second, because of limited computing resources, our study precludes a full scale factorial experiment which, if otherwise taken, would have provided significant insights into the effects of different p a r a m e t e r combinations. For the same reason, the number of facilities is limited to 30 in our experiments. To demonstrate the effectiveness of G A in laying out a large number of facilities, the problem size must be increased. Third, the use of numerical clustering is shown to be relatively efficient but is surely not the only way to construct a slicing tree structure. Other alternatives such as min cut and Steiner tree should be considered. Fourth, the geometry of a facility can be extended beyond rectangle to other shapes. Last, the procedure can be extended to support an interactive dialogue with the designer. For instance, the system may maintain a list of good solutions encountered so far. The designer could halt the search at any time and inspect the entries in the list. He may also 'freeze' some cut operators in the sequence if he thinks some parts of a layout (i.e., subtrees) are good enough and need not be changed. For example, each column of layouts shown in Figure 28 is obtained by fixing the first operator. Similarly, each row in the diagram represent layouts with their second operator fixed. In addition, the designer may decide to reduce the crossover and mutation rate if existing solutions are acceptable and a fast convergence is needed. As a result, schemas with high fitness will have a higher sampling rate. The overall effect is to reproduce an increasing number of representatives of the best schemas in each generation with the extreme case that all solutions of a population are identical.

References Anderbcrg, M.R. (1973), Cluster Analysis for Applications, Academic Press, New York. Armour, G.C., and Buffa, E.S. (1963), "A heuristic algorithm and simulation approach to the relative location of facilities", Management Science 9, 294-309. Brooker, L.B., Goldberg, D.E., and Holland, J.H. (1989), "Classifier systems and genetic algorithms", Artificial httelligence 40, 235 282. Chu, C.H. (1989), "Clustering analysis in manufacturing cellular formation", OMEGA 17, 289-295.

346

K.Y. Tam / Genetic algorithms, optim&ation and facility layout design

DeJong, K.A. (1975), "Analysis of the behavior of a class of genetic adaptive systems", Ph.D. Dissertation, Department of Computer and Communication Sciences, University of Michigan, Ann Arbor, MI. DeJong, K.A. (1988), "Learning with genetic algorithms: An overview", Machine Learning 3, 121-138. Goldberg, D.E. (1989), Genetic Algorithms in Search, Optimization, and Machine Learning, Addison-Wesley, Reading, MA. Grefenstette, J.J. (1986), "Optimization of control parameters for genetic algorithms", IEEE Transactions on Systems, Man, and Cybernetics 16, 122-128. Grefenstette, J.J., and Baker, J.E. (1989), "How genetic algorithms: A critical look at implicit parallelism", in: Proceedings of the Third International Conference on Genetic Algorithms, Morgan Kaufmann, Los Altos, CA, 20-27. Heragu, S.S., and Kusiak, A. (1988), "Machine layout problem in flexible manufacturing systems", Operations Research 36, 258-268. Holland, J.H. (1973), "Genetic algorithms and the optimal allocations of trials", SlAM Journal of Computing 2, 88-105. Holland, J.H. (1975), Adaptation in Natural and Artificial Systems, The University of Michigan Press, Ann Arbor, MI. Holland, J.H. (1980), "Adaptive algorithms for discovering and using general patterns in growing knowledge-bases", International Journal of Policy Analysis and Information Systems 4, 217-240. Kirkpatrick, C., Gelatt, C., and Vecchi, M. (1983), "Optimization by simulated annealing", Science 220, 671-680. Laarhoven, P.J.M. and Aarts, E.H.L. (1987), Simulated Annealing." Theory and Applications, Kluwer, Boston, MA. Lee, R.C., and Moore, J.M. (1967), "CORELAP - Computerized Relationship Layout Planning", Industrial Engineering 18, 195-200. Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A., and Teller, E. (1953), "Equation of state calculations by fast computing machines", Journal of Chemical Physics 21, 1087-1092. Moore, J.M. (1974), "Computer aided facilities design: An international survey", International Journal of Production Research 12, 21-44. Nugent, C.E,, Vollman, T.E., and Ruml, J. (1968), "An experimental comparison of techniques for the assignment of facilities to locations", Operations Research 16, 150-173. Otten, R.H.J.M., and Ginneken, L.P.P.P. (1989), The Annealing Algorithm, Kluwer, Boston, MA. Sahni, S., and Gonzalez, T. (1976), "P-complete approximation problem", Journal of ACM 23, 555-565. Schaffer, J.D., Caruana, R.A., Eshelman, L.J., and Das, R. (1989), "A study of control parameters affecting online performance of genetic algorithms for function optimization", in: Proceedings of the Third International Conference on Genetic Algorithms, Morgan Kaufmann, Los Altos, CA, 51-60. Seehof, J.M., and Evans, W.O. (1967), "Automated layout design program", Industrial Engineering 18, 690-695. Skiscim, C., and Golden, B. (1986), "Using simulated annealing to solve routing and location problems", Naval Research Logistics Quarterly 33, 261-273. Tam, K.Y. (1991), "A simulated annealing algorithm for allocating space to manufacturing cells", International Journal of Production Research (forthcoming). Tompkins, J.A., and White, J.A. (1984), Facilities Planning, Wiley, New York. Ward, J.H. (1963), "Hierarchical grouping to optimize an objective function", Journal of the American Statistical Association 58, 236-244.