Numerical study of variations of airflow induced by a moving automatic guided vehicle in a cleanroom

(1)

This article was downloaded by: [National Chiao Tung University 國立交通大學] On: 27 April 2014, At: 22:21

Publisher: Taylor & Francis

Informa Ltd Registered in England and Wales Registered Number: 1072954 Registered office: Mortimer House, 37-41 Mortimer Street, London W1T 3JH, UK

Journal of the Chinese Institute of Engineers

Publication details, including instructions for authors and subscription information: http://www.tandfonline.com/loi/tcie20

Numerical study of variations of airflow induced by a

moving automatic guided vehicle in a cleanroom

Suh‐Jenq Yang a_{, Shih‐Fa Chen}b_{& Wu‐Shung Fu}b a

Department of Industrial Engineering and Management , Nan Kai College , Nantou, Taiwan 542, R.O.C.

b

Department of Mechanical Engineering , National Chiao Tung University , Hsinchu, Taiwan 300, R.O.C.

Published online: 03 Mar 2011.

To cite this article: Suh‐Jenq Yang , Shih‐Fa Chen & Wu‐Shung Fu (2002) Numerical study of variations of airflow induced

by a moving automatic guided vehicle in a cleanroom, Journal of the Chinese Institute of Engineers, 25:1, 67-75, DOI: 10.1080/02533839.2002.9670681

To link to this article: http://dx.doi.org/10.1080/02533839.2002.9670681

PLEASE SCROLL DOWN FOR ARTICLE

Taylor & Francis makes every effort to ensure the accuracy of all the information (the “Content”) contained in the publications on our platform. However, Taylor & Francis, our agents, and our licensors make no

representations or warranties whatsoever as to the accuracy, completeness, or suitability for any purpose of the Content. Any opinions and views expressed in this publication are the opinions and views of the authors, and are not the views of or endorsed by Taylor & Francis. The accuracy of the Content should not be relied upon and should be independently verified with primary sources of information. Taylor and Francis shall not be liable for any losses, actions, claims, proceedings, demands, costs, expenses, damages, and other liabilities whatsoever or howsoever caused arising directly or indirectly in connection with, in relation to or arising out of the use of the Content.

This article may be used for research, teaching, and private study purposes. Any substantial or systematic reproduction, redistribution, reselling, loan, sub-licensing, systematic supply, or distribution in any

form to anyone is expressly forbidden. Terms & Conditions of access and use can be found at http:// www.tandfonline.com/page/terms-and-conditions

(2)

NUMERICAL STUDY OF VARIATIONS OF AIRFLOW INDUCED

BY A MOVING AUTOMATIC GUIDED VEHICLE IN A CLEANROOM

Suh-Jenq Yang

Department of Industrial Engineering and Management Nan Kai College

Nantou, Taiwan 542, R.O.C. Shih-Fa Chenand Wu-Shung Fu*

Department of Mechanical Engineering National Chiao Tung University

Hsinchu, Taiwan 300, R.O.C.

Key Words: cleanroom, AGV, recirculation zones, FEM.

ABSTRACT

The variations of airflow induced by a moving automatic guided vehicle (AGV) in a vertical laminar flow cleanroom are studied numerically. From a viewpoint of fluid mechanics, the characteristic of the variations of the airflow induced by a moving object is dynamic and is classified as a moving boundary problem. A Galerkin finite element formulation with an arbitrary Lagrangian-Eulerian (ALE) ki-nematic description method is adopted to analyze this problem. Three different moving velocities of the AGV under Reynolds number Re=500 and two different positions of the wafer cassette are considered. The results show that the formation of recirculation zones, which are dis-advantageous for removing contaminants, is remarkably dependent on the moving velocity of the AGV and the position of the wafer cassette.

*Correspondence addressee

I. INTRODUCTION

Recently, accompanying state-of-the-art semi-conductor techniques, wafers have increased in both size and weight. It has become difficult for an op-erator without any auxiliary tools to carry the larger and heavier wafers. Thus, an effective transport tool is required to save labor and decrease wafer damage and contamination. An AGV is commonly used in a cleanroom for transporting wafers.

In the past, most studies about AGVs were fo-cused on the characteristics, system control, and path planning. Among studies devoted to the characteris-tics of the AGV, McClelland (1986) described the potential advantages of the utilization of the AGV in

the semiconductor industry to prevent contamination. Lee et al. (1996) used five heuristic rules to deal with the load selection of the AGV. Rajotia et al. (1998) employed analytical and simulation modeling to de-termine the number of AGVs needed. Referring to the literature of process control, Zaremba et al. (1997) adopted max-algebra formalism for the robust distrib-uted control of the AGV. Funabiki and Kodera (1998) proposed a steering control strategy of the AGV with successive learning neural networks using a real-time tuning function. Kim et al. (1999) utilized a dead-lock prevention algorithm to describe the dispatch-ing problem based on workload balancdispatch-ing for the control of the AGV. Concerning research on path planning, Lin and Wang (1997) proposed a fuzzy

(3)

68 Journal of the Chinese Institute of Engineers, Vol. 25, No. 1 (2002)

approach, which included sensor modeling and trap recovering to prevent collision with unknown mov-ing obstacles. Oboth et al. (1999) adopted a large-scale simulation of a dynamic batch type to present a route-generation technique that provided conflict-free routes for multiple AGVs. Sekine et al. (1999) used the imbedded Markov chain method to deal with the analysis of traffic congestion among AGVs.

As for the other important issue, the effect of the AGV on the motions of the airflow in the cleanroom has been relatively less studied. Kanayama et al. (1993) adopted a numerical method to analyze a t w o - d i m e n s i o n a l f l o w a r o u n d a n A G V i n a cleanroom. The results showed that the wafer might not be contaminated by particles from the floor. However, to facilitate the analysis, the study men-tioned above regarded the moving AGV as a station-ary one, which resulted in the phenomena of the air-flow being rather different from what actually occurs. Consequently, the aim of this study is to numeri-cally simulate the motions of the airflow affected by a moving AGV in the cleanroom. Because of the movement of the AGV, this subject is classified as a moving boundary problem, which is difficult to ana-lyze solely using either the Lagrangian or Eulerian kinematic description method. An appropriate kine-matic description method of the ALE method, which combines the characteristics of the Lagrangian and Eulerian kinematic description methods, is adopted to describe this problem. In the ALE method, the computational meshes may move with the fluid (Lagrangian), be held fixed (Eulerian), or be moved in any other prescribed way. The details of the kine-matic theory of the ALE method are given in Hughes

et al.. (1981), Donea et al. (1982), and Ramaswamy

(1990). A Galerkin finite element method with mov-ing meshes and an implicit difference scheme, deal-ing with the time terms, is used to solve the govern-ing equations. Three different movgovern-ing velocities of the AGV and two different positions of the wafer cas-sette are considered. The results show that the for-mation of recirculation zones around the AGV and wafer cassette is remarkably dependent on the mov-ing velocity of the AGV and the position of the wafer cassette. These phenomena are quite different from those shown in Kanayama et al. (1993) in which the moving AGV was regarded as a stationary one.

II. PHYSICAL MODEL

A model of a two-dimensional vertical laminar flow cleanroom sketched in Fig. 1 is used. The height and width of the cleanroom are h(=h1+h2+h3) and w,

respectively. Two workbenches, with height h2 and

width w1, are set on the right and left sides of the

cleanroom, respectively. An AGV, with height h5 and

width w8(=w3+w4+w5), moves between the two

workbenches. The wafer cassette, with height h4 and

width w4, is put on the top surface of the AGV. The

inlet air, velocity v0, flows into the cleanroom at the

inlet section AB. Initially (t=0), the distance from the left workbench to the AGV is w2, and the AGV is

stationary and the airflow flows steadily. As the time

t>0, the AGV moves toward the right workbench with

a constant velocity ub. The interaction between the airflow and AGV affects the behavior of the airflow in the cleanroom. This subject is classified as a mov-ing boundary problem and becomes a time-dependent problem. Thus, the ALE method is properly utilized to analyze this problem.

In order to facilitate the analysis, the following assumptions are made.

(1) The flow field is two-dimensional, incompressible and laminar.

(2) The properties of the fluid are constant and the effect of gravity is neglected.

(3) The no-slip condition is held on the interfaces between the airflow and AGV.

According to the characteristic scales of w1, v0

and ρv02, the dimensionless variables are defined as

follows: X = x_w 1, Y = y w1, U = uv0, V = v_v 0, U = uv0, Ub= ub v0 , Re =v0w1 ν , τ=tv0 w1 , P = p – p_∞ ρv02 , (1)

Fig. 1 Physical model.

(4)

where u is the mesh velocity in the x-direction. Based upon the above assumptions and dimen-sionless variables, the dimendimen-sionless ALE governing equations are expressed as the following equations: continuity equation ∂U ∂X + ∂V ∂Y = 0 (2) momentum equations ∂U ∂τ + (U – U)∂∂UX + V∂ U ∂Y = –∂ P ∂X + 1Re( ∂2 U ∂X2 + ∂2 U ∂Y2) (3) ∂V ∂τ + (U – U)∂∂VX+ V∂ V ∂Y = –∂ P ∂Y + 1Re( ∂2 V ∂X2+ ∂2 V ∂Y2) (4)

As the time τ>0, the boundary conditions are as follows:

on the wall surfaces of the cleanroom

U=V=0, (5)

on the airflow inlet section AB

U=0, V=-1.0, (6)

on the airflow outlet section EF

∂U

∂Y =∂ V

∂Y = 0 , (7)

on the interfaces of the AGV and the airflow

U=Ub, V=0. (8)

III. NUMERICAL METHOD

A Galerkin finite element formulation with mov-ing meshes and an implicit scheme, dealmov-ing with the time terms, is adopted to solve the governing Eqs. (2)-(4). The Newton-Raphson method is used to lin-e a r i z lin-e t h lin-e n o n l i n lin-e a r t lin-e r m s i n t h lin-e m o m lin-e n t u m equations, and the pressure is eliminated from the momentum equations using the penalty function model (Reddy and Gartling, 1994). Nine node quad-rilateral elements are utilized to discretize the prob-lem domain and the velocity terms are approximated using the quadratic shape functions. The discretiza-tion processes of the governing equadiscretiza-tions are similar to the one used in Fu et al. (1990). Then, the mo-mentum Eqs. (3)-(4), can be expressed as follows:

([A](e)+ [K](e)+λ[L](e)){q}_τ+∆τ (e)

Σ

1 n e =

Σ

{f}(e) 1 n e , (9) where ({q}_τ(e)₊_∆τ)T= U₁, U₂, , U₉, V₁, V₂, , V₉ τ+∆τ m + 1 , (10) [A](e)

consists of the (m)th iteration values of U and V at the time τ+∆τ,

[K](e)

consists of the U and time differential terms,

[L](e)

consists of the penalty function terms, {f}(e)

consists of the known values of U and V at the time τ and (m)th iteration values of

U and V at the time τ+∆τ.

In Eq. (9), the terms with the penalty parameter λ are integrated by 2×2 Gaussian quadrature, and the o t h e r t e r m s a r e i n t e g r a t e d b y 3×3 G a u s s i a n quadrature. The value of penalty parameter used in this study is 106

, and the frontal method solver Tay-lor and Hughes (1981) is used to solve Eq. (9).

The mesh velocity U is assumed to be linearly distributed and inversely proportional to the distance between the computational meshes and the AGV. To prevent the computational nodes in the vicinity of the AGV slipping away from the boundary layer, the mesh velocities adjacent to the AGV are expediently as-signed equal to the velocity of the AGV.

A brief outline of the solution procedure follows: (1) Determine the optimal mesh distribution and

num-ber of elements and nodes.

(2) Solve the values of U and V at the steady state and regard them as the initial values.

(3) Determine the time increment ∆τ and the mesh ve-locities U at every node.

(4) Update the coordinates of the nodes and examine the determinant of the Jacobian transformation matrix to ensure the one-to-one mapping is satis-fied during the Gaussian quadrature numerical i n t e g r a t i o n ; o t h e r w i s e , e x e c u t e t h e m e s h reconstruction.

(5) Solve Eq. (9) until the following criteria for con-vergence are satisfied:

φm + 1 –φm φm + 1 τ+∆τ < 10– 3, where φ=U, V (11)

(6) Continue the next time step calculation until the assigned position of the AGV is reached.

IV. RESULTS AND DISCUSSION

Four cases with different moving velocities and the dimensionless geometric parameters are tabulated in Table 1. Cases 1, 2, and 3, which corre-spond to the moving velocities Ub=0.5, 2.0, and 4.0,

(5)

respectively, have the same geometric condition. Case 4 has the same moving velocity as case 2, but the position of the wafer cassette in case 4 is on the left top surface of the AGV. The variations of the airflow induced by the movement of the AGV under

Re=500 are examined in detail.

For matching the boundary conditions at the in-let and outin-let of the cleanroom mentioned earlier, the lengths from the inlet and outlet to the AGV are de-termined by numerical tests and are equal to 8.0 and 36.0, respectively. To obtain an optimal computa-tional mesh, three different nonuniform distributed elements 3116, 3456, and 3892 (corresponding to 12572, 14124, and 15892 nodes, respectively) are

tested for case 1 as the AGV is stationary. The re-sults of the distributions of velocities U and V along the lines MM′ and NN′, as indicated in Fig. 1, are shown in Fig. 2. According to the results of the mesh tests, the computational mesh with 3892 elements is adopted for cases 1, 2, and 3. Similarly, the nonuni-form distribution of 4000 elements corresponding to 16332 nodes is adopted for case 4. Moreover, an im-plicit scheme is employed to deal with the time dif-ferential terms of the governing equations. The time step ∆τ=0.005 is chosen for all cases in this study.

To illustrate the variations of the airflow in more detail, we only presented the velocity vectors around the AGV and wafer cassette. Besides, the

Table 1 The dimensionless geometric parameters of the cleanroom for four different cases

Ub H1 H2 H3 H4 H5 H6 W W1 W2 W3 W4 W5 W6

case 1 0.5 7 4 34 1 1 36 22 1 4 1 2 1 12

case 2 2.0 7 4 34 1 1 36 22 1 4 1 2 1 12

case 3 4.0 7 4 34 1 1 36 22 1 4 1 2 1 12

case 4 2.0 7 4 34 1 1 36 22 1 4 0.25 2 1.75 12

Fig. 2 Comparison of the distributions of the velocities U and V along the lines of MM′ and NN′ at the steady state and Re=500 for various meshes (a) X-U, (b) X-V, (c) Y-U, (d) Y-V

(6)

velocity vectors shown in the following figures are scaled relative to the maximum velocity in the flow field.

Figure 3 shows the transient developments of the velocity vectors around the AGV and wafer cassette for case 1. In this case, the moving velocity

Ub of the AGV is 0.5. At the time τ=0.0, as shown in Fig. 3(a), the AGV is stationary and the distance w2

is equal to 4.0. The airflow flows downward, steadily. A large recirculation zone is observed near the bot-tom surface of the AGV. As the time τ>0.0, the AGV moves toward the right workbench with a constant

moving velocity Ub=0.5. Thus, the variations of the flow field become the transient state. Because the velocity of the inlet airflow is larger than the moving velocity of the AGV, the AGV seems to be an ob-stacle in the way of the airflow. At the time τ=4.0 (Fig. 3(b)), the AGV presses the airflow near the front surface of the AGV, and the direction of the airflow is forced to change and turn into the bottom surface of the AGV. As for the airflow near the rear region of the AGV, the inlet airflow simultaneously replen-ishes the vacant space induced by the movement of the AGV. Consequently, new recirculation zones are formed around the rear region of the AGV and cassette. Particles might suspend in these recircula-tion zones near the AGV and cassette and are diffi-cult for the airflow to remove. As the time increases, as shown in Figs. 3(c)-(d), the recirculation zones around the rear region of the AGV and cassette en-large gradually, and the variations of the airflow around the AGV become more complex. These phenomena are disadvantageous to removal of par-ticles from the cleanroom. Furthermore, because of the moving of the AGV, some of the airflow flow-ing through the rear region of the AGV flows toward the bottom surface of the AGV. As a result, the air-flow might sweep particles from the floor to the work areas and this is disadvantageous for the semi-conductor manufacturing process. The results are re-markably different from those shown in Kanayama

et al. (1993), which regarded the moving AGV as a

stationary one.

Figure 4 shows the transient developments of the velocity vectors around the AGV and wafer cas-sette for case 2. In this case, the AGV moves toward the right workbench with a constant velocity Ub= 2.0, which is much faster than that of the above case. In the first stage of the transient state, as shown in Fig. 4(b), the variations of the airflow are similar to the above case. As the time increases, as shown in Figs. 4(c)-(d), because the moving velocity of the AGV is faster than that of case 1, the recirculation zones around the AGV and cassette are larger in this case.

The transient developments of the velocity vec-tors around the AGV and wafer cassette for case 3, in which the moving velocity of the AGV is equal to 4.0, are shown in Fig. 5. Basically, the variations of the flow field are similar to case 2. However, the strengths of recirculation zones around the rear re-gion of the AGV and cassette are stronger than those of the above cases; that is disadvantageous for re-moving particles. Besides, because of the friction between the AGV and rail, the increase of the mov-ing velocity of the AGV will produce more particles, which results in serious contamination on the products.

Fig. 3 The transient developments of the velocity vectors around the AGV and wafer cassette for case 1 (Ub=0.5) (a) τ=0.0,

(b) τ=4.0, (c) τ=8.0, (d) τ=16.0

(7)

As shown in Fig. 6, there are the transient de-velopments of the velocity vectors around the AGV and wafer cassette for case 4, in which the moving velocity of the AGV is the same as case 2, but the wafer cassette is near the left side of the top surface of the AGV. At the time τ=1.0, the recirculation zones around the rear region of the AGV and cassette are similar to case 2; however, the recirculation zone near the cassette is larger than that of case 2. At the time τ=2.0, the recirculation zones around the rear region of the AGV and cassette enlarge gradually. As the time increases, the recirculation zones around the rear

(b) τ=0.5, (c) τ=1.0, (d) τ=2.0

region of the AGV and cassette shrink gradually; these phenomena are somewhat different from those in case 2. Therefore, the position of the wafer cassette near the lateral side seems to be better than that in the cen-tral region of the top surface of the AGV.

V. CONCLUSIONS

The variations of the airflow induced by the movement of an AGV in a cleanroom are investigated numerically. The main conclusions can be summa-rized as follows:

(b) τ=1.0, (c) τ=2.0, (d) τ=4.0

(8)

1. The recirculation zones, in which particles are easily trapped, are observed around the AGV and wafer cassette as the AGV moves. The results are different from studies regarding a moving AGV as a stationary one.

2. The lower the moving velocity of the AGV, the smaller the recirculation zones formed.

3. From a viewpoint of contamination control, the po-sition of the wafer cassette near either side of the top surface of the AGV is better than that in the central region of the top surface of the AGV.

ACKNOWLEDGEMENT

The support of this work by the National Sci-ence Council of Taiwan, R.O.C., under contract NSC90-2626-E-252-003 is gratefully acknowledged.

NOMENCLATURE

h the height of the cleanroom [m]

H dimensionless height of the cleanroom

h1 the distance from the workbench to the ceiling

of the cleanroom [m]

H1 dimensionless distance from the workbench to

the ceiling of the cleanroom

h2 the height of the workbench [m]

H2 dimensionless height of the workbench

h3 the distance from the workbench to the outlet

[m]

H3 dimensionless distance from the workbench to

the outlet

h4 the height of the cassette [m]

H4 dimensionless height of the cassette

h5 the height of the AGV [m]

H5 dimensionless height of the AGV

h6 the distance from the AGV to the outlet [m]

H6 dimensionless distance from the AGV to the

outlet p pressure [N/m2 ] p∞ reference pressure [N/m2 ] P dimensionless pressure Re Reynolds number t time [s]

u, v velocities of the airflow in x and y directions [m/s]

U,V dimensionless velocities of the airflow in X and

Y directions

ub the moving velocity of the AGV [m/s]

Ub dimensionless moving velocity of the AGV

u the mesh velocity in x-direction [m/s]

U dimensionless mesh velocity in X-direction

v0 the airflow inlet velocity [m/s]

V0 dimensionless airflow inlet velocity

w the width of the cleanroom [m]

W dimensionless width of the cleanroom

w1 the width the of workbench [m]

W1 dimensionless width of the workbench

w2 the distance from the left workbench to the

AGV [m]

W2 dimensionless distance from the left

work-bench to the AGV

w4 the width of the cassette [m]

W4 dimensionless width of the cassette

w6 the distance from the AGV to the right

work-bench [m]

W6 dimensionless distance from the AGV to the

right workbench Fig. 6 The transient developments of the velocity vectors around

the AGV and wafer cassette for case 4 (Ub=2.0) (a) τ=0.0,

(b) τ=1.0, (c) τ=2.0, (d) τ=4.0

(9)

x, y Cartesian coordinates [m]

X, Y dimensionless Cartesian coordinates

Greek λ penalty parameter ν kinematic viscosity [m2 /s] ρ density [kg/m3 ] τ dimensionless time Other [] matrix {} column vector <> row vector | | absolute value REFERENCES

1. Donea, J., Giuliani, S., and Halleux, J. P., 1982, “An Arbitrary Lagrangian-Eulerian Finite Element Method for Transient Dynamic Fluid-structure Interactions,” Computer Methods in

Ap-plied Mechanics and Engineering, Vol. 33, pp.

689-723.

2. F u , W . S . , K a u , T . M . , a n d S h i e h , W . J . , 1990, “Transient Laminar Natural Convection in an Enclosure from Steady Flow State to Station-ary State,” Numerical Heat Transfer, Part A, Vol. 18, pp. 189-211.

3. Funabiki, S. and Kodera, K., 1998, “A Steering Control of Automated Guided Vehicles by the Neural Networks Using a Real-time Function,”

Transactions of the Institute of Electrical Engi-neers Japan, Part D, Vol. 118, No. 5, pp.

605-610.

4. Hughes, T. J. R., Liu, W. K., and Zimmermann, T. K., 1981, “Lagrangian-Eulerian Finite Element Formulation for Incompressible Viscous Flows,”

Computer Methods in Applied Mechanics and Engineering, Vol. 29, pp. 329-349.

5. Kanayama, H., Toshigami, K., Tashiro, Y., Tabata, M., and Fujima, S., 1993, “Finite Element Analysis of Air Flow Around an Automatic Guided Vehicle,” Journal of Wind Engineering

and Industrial Aerodynamics, Vol. 46-47, pp.

801-810.

6. Kim, C. W., Tanchoco, J. M. A., and Koo, P. H.,

1999, “AGV Dispatching Based on Workload Balancing,” International Journal of Production

Research, Vol. 37, No. 17, pp. 4053-4066.

7. Lee, J., Tangjarukij, M., and Zhu, Z., 1996, “Load Selection of Automated Guided Vehicles in Flex-ible Manufacturing Systems,” International

Jour-nal of Production Research, Vol. 34, No. 12, pp.

3383-3400.

8. Lin, C. H. and Wang L. L., 1997, “Intelligent Collision Avoidance by Fuzzy Logic Control,”

Robotics and Autonomous Systems, Vol. 20, No.

1, pp. 61-83.

9. McClelland, S., 1986, “The Cleanest Robots in the World’ (Clean Room Robots),” The Industrial

Robot, Vol. 13, No. 4, pp. 217-220.

10. Oboth, C., Batta, R., and Karwan M., 1999, “Dy-namic Conflict-Free Routing of Automated Guided Vehicles,” International Journal of

Pro-duction Research, Vol. 37, No. 9, pp. 2003-2030.

11. Rajotia, S., Shanker, K., and Batra, J. L., 1998, “Determination of Optimal AGV Fleet Size for an FMS,” International Journal of

Produc-tion Research, Vol. 36, No. 5, pp. 1177-1198.

12. Ramaswamy, B., 1990, “Numerical Simulation of Unsteady Viscous Free Surface Flow,” Journal

of Computational Physics, Vol. 90, pp. 396-430.

13. Reddy, J. N. and Gartling, D. K., 1994, The Fi-nite Element Method in Heat Transfer and Fluid Dynamics, CRC Press. Inc., Ann. Arbor, chapter 4.

14. Sekine, E., Hamamatsu, Y., and Kongouji, T., 1999, “Analysis of Traffic Congestion in AGV System,” Transactions of the Institute of

Electri-cal Engineers Japan, Part D, Vol. 119, No. 4, pp.

515-522.

15. Taylor, C. and Hughes, T. G., 1981, Finite

Ele-m e n t P r o g r a Ele-m Ele-m i n g o f t h e N a v i e r - S t o k e s Equations, Pineridge Press Ltd., U.K., chapter 6.

16. Zaremba, M. B., Obuchowicz, A., Banaszak, Z. A., and Jedrzejek, K. J., 1997, “A Max-Alge-bra Approach to the Robust Distributed Control of Repetitive AGV Systems,” International

Jour-nal of Production Research, Vol. 35, No. 10, pp.

2667-2687.

Manuscript Received: Apr. 09, 2001 Revision Received: Aug. 07, 2001 and Accepted: Oct. 02, 2001

(10)

!"#$%&'()*+,-./,01

!"#

!"#$%$&'()* !"#$ !"#$%&'#( !"#$%&'()*+,-./012/34567890:; !"#$%&'()*+,-./,0123456789 !" !"# !"#$%&'()*+,-./0123456789 !"#$%&'()*+ Re=500 !"#$%&'()*+,- !"#$%&'()*+,-.,/01234567 8!"#$%&'()*+,- 9: !"#$%&'()*+,-#./-0123456789:;9< !"#$%&'()*+ !",-./012%34*5678%9 !"#

!"#$%&'()*+%,-.%/01234