ISSN 1854-6250 APEM journal Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 ü Published by PEI apem-journal.org University of Mari bor Advances in Production Engineering & Management Identification Statement APEM University of Maribor ISSN 1854-6250 | Abbreviated key title: Adv produc engineer manag | Start year: 2006 ISSN 1855-6531 (on-line) Published quarterly by Production Engineering Institute (PEI), University of Maribor Smetanova ulica 17, SI - 2000 Maribor, Slovenia, European Union (EU) Phone: 00386 2 2207522, Fax: 00386 2 2207990 Language of text: English APEM homepage: apem-journal.org University homepage: www.um.si APEM Editorial Editor-in-Chief Miran Brezocnik editor@apem-journal.org, info@apem-journal.org University of Maribor, Faculty of Mechanical Engineering Smetanova ulica 17, SI - 2000 Maribor, Slovenia, EU Desk Editor Tomaz irgolic deski@apem-journal.org Website Technical Editor Lucija Brezocnik lucija.brezocnik@um.si Editorial Board Members Eberhard Abele, Technical University of Darmstadt, Germany Bojan Acko, University of Maribor, Slovenia Joze Balic, University of Maribor, Slovenia Agostino Bruzzone, University of Genoa, Italy Borut Buchmeister, University of Maribor, Slovenia Ludwig Cardon, Ghent University, Belgium Nirupam Chakraborti, Indian Institute of Technology, Kharagpur, India Edward Chlebus, Wroclaw University of Technology, Poland Franci Cus, University of Maribor, Slovenia Igor Drstvensek, University of Maribor, Slovenia Illes Dudas, University of Miskolc, Hungary Mirko Ficko, University of Maribor, Slovenia Vlatka Hlupic, University of Westminster, UK David Hui, University of New Orleans, USA Pramod K. Jain, Indian Institute of Technology Roorkee, India Isak Karabegovic, University of Bihac, Bosnia and Herzegovina Janez Kopac, University of Ljubljana, Slovenia Iztok Palcic, University of Maribor, Slovenia Krsto Pandza, University of Leeds, UK Andrej Polajnar, University of Maribor, Slovenia Antonio Pouzada, University of Minho, Portugal Rajiv Kumar Sharma, National Institute of Technology, India Katica Simunovic, J. J. Strossmayer University of Osijek, Croatia Daizhong Su, Nottingham Trent University, UK Soemon Takakuwa, Nagoya University, Japan Nikos Tsourveloudis, Technical University of Crete, Greece Tomo Udiljak, University of Zagreb, Croatia Ivica Veza, University of Split, Croatia Limited Permission to Photocopy: Permission is granted to photocopy portions of this publication for personal use and for the use of clients and students as allowed by national copyright laws. This permission does not extend to other types of reproduction nor to copying for incorporation into commercial advertising or any other profit-making purpose. Subscription Rate: 120 EUR for 4 issues (worldwide postage included); 30 EUR for single copies (plus 10 EUR for postage); for details about payment please contact: info@apem-journal.org Cover and interior design: Miran Brezocnik Printed: Tiskarna Kostomaj, Celje, Slovenia Subsidizer: The journal is subsidized by Slovenian Research Agency Statements and opinions expressed in the articles and communications are those of the individual contributors and not necessarily those of the editors or the publisher. No responsibility is accepted for the accuracy of information contained in the text, illustrations or advertisements. Production Engineering Institute assumes no responsibility or liability for any damage or injury to persons or property arising from the use of any materials, instructions, methods or ideas contained herein. Copyright © 2016 PEI, University of Maribor. All rights reserved. Advances in Production Engineering & Management is indexed and abstracted in the WEB OF SCIENCE (maintained by THOMSON REUTERS): Science Citation Index Expanded, Journal Citation Reports - Science Edition, Current Contents - Engineering, Computing and Technology • Scopus (maintained by Elsevier) • Inspec • EBSCO: Academic Search Alumni Edition, Academic Search Complete, Academic Search Elite, Academic Search Premier, Engineering Source, Sales & Marketing Source, TOC Premier • ProQuest: CSA Engineering Research Database -Cambridge Scientific Abstracts, Materials Business File, Materials Research Database, Mechanical & Transportation Engineering Abstracts, ProQuest SciTech Collection • TEMA (DOMA) • The journal is listed in Ulrich's Periodicals Directory and Cabell's Directory University of Maribor Production Engineering Institute (PEI) Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 | pp 255-380 Contents Scope and topics 258 A combined zone-LP and simulated annealing algorithm for unequal-area 259 facility layout problem Xiao, Y.J.; Zheng, Y.; Zhang, L.M.; Kuo, Y.H. A new multi-objective Jaya algorithm for optimization of modern machining processes 271 Rao, R.V.; Rai, D.P.; Ramkumar, J.; Balic, J. Finite element method for optimum design selection of carport structures under 287 multiple load cases Ozkal, F.M.; Cakir, F.; Arkun, A.K. Applying multi-phase particle swarm optimization to solve bulk cargo port 299 scheduling problem Tang, M.; Gong, D.; Liu, S.; Zhang, H. A case-based reasoning approach for non-traditional machining processes selection 311 Boral, S.; Chakraborty, S. Dimensional accuracy of camera casing models 3D printed on Mcor IRIS: A case study 324 Mandic, M.; Galeta, T.; Raos, P.; Jugovic, V. Effect of delayed differentiation on a multiproduct vendor-buyer integrated 333 inventory system with rework Chiu, Y.-S.P.; Kuo, J.-S.; Chiu, S.W.; Hsieh, Y.-T. Experimental modeling of fluid pressure during hydroforming of welded plates 345 Karabegovic, E.; Poljak, J. An integrated lean approach to Process Failure Mode and Effect Analysis (PFMEA): 355 A case study from automotive industry Banduka, N.; Veza, I.; Bilic, B. Multi-objective optimization of the turning process using a Gravitational Search Algorithm 366 (GSA) and NSGA-II approach Klancnik, S.; Hrelja, M.; Balic, J.; Brezocnik, M. Calendar of events 377 Notes for contributors 379 Journal homepage: apem-journal.org ISSN 1854-6250 (print) ISSN 1855-6531 (on-line) ©2016 PEI, University of Maribor. All rights reserved. 257 Scope and topics Advances in Production Engineering & Management (APEM journal) is an interdisciplinary refer-eed international academic journal published quarterly by the Production Engineering Institute at the University of Maribor. The main goal of the APEM journal is to present original, high quality, theoretical and application-oriented research developments in all areas of production engineering and production management to a broad audience of academics and practitioners. In order to bridge the gap between theory and practice, applications based on advanced theory and case studies are particularly welcome. For theoretical papers, their originality and research contributions are the main factors in the evaluation process. General approaches, formalisms, algorithms or techniques should be illustrated with significant applications that demonstrate their applicability to real-world problems. Although the APEM journal main goal is to publish original research papers, review articles and professional papers are occasionally published. Fields of interest include, but are not limited to: Additive Manufacturing Processes Advanced Production Technologies Artificial Intelligence in Production Assembly Systems Automation Computer-Integrated Manufacturing Cutting and Forming Processes Decision Support Systems Discrete Systems and Methodology e-Manufacturing Evolutionary Computation in Production Fuzzy Systems Human Factor Engineering, Ergonomics Industrial Engineering Industrial Processes Industrial Robotics Intelligent Manufacturing Systems Joining Processes Knowledge Management Logistics Machine Learning in Production Machine Tools Machining Systems Manufacturing Systems Materials Science, Multidisciplinary Mechanical Engineering Mechatronics Metrology Modelling and Simulation Numerical Techniques Operations Research Operations Planning, Scheduling and Control Optimisation Techniques Project Management Quality Management Risk and Uncertainty Self-Organizing Systems Statistical Methods Supply Chain Management Virtual Reality in Production 258 Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 | pp 259-270 http://dx.doi.Org/10.14743/apem2016.4.225 ISSN 1854-6250 Journal home: apem-journal.org Original scientific paper A combined zone-LP and simulated annealing algorithm for unequal-area facility layout problem Xiao, Y.J.ab, Zheng, Y.c, Zhang, L.M.d*, Kuo, Y.H.e aJiangsu Key Laboratory of Modern Logistics, School of Marketing and Logistic Management, Nanjing University of Finance and Economics, Nanjing , China bBusiness school, Nanjing University, Nanjing, China cCollege of Automobile and Traffic Engineering, Nanjing Forestry University, Nanjing, China dSchool of Management and Engineering, Nanjing University, Nanjing, China eStanley Ho Big Data Decision Analytics Research Centre, The Chinese University of Hong Kong, Shatin, New Territories, Hong Kong A B S T R A C T A R T I C L E I N F O Facility layout problem (FLP) is one of well-known NP-hard problems and has been demonstrated to be useful in enhancing the productivity of manufacturing systems in practice. This paper focuses on the unequal-area FLP (UA-FLP) whose goal is to locate departments with different areas within a given facility so as to minimize the total material handling cost. A novel approach, which we call a combined zone-linear programming (zone-LP) and simulated annealing algorithm, is developed for solving the UA-FLP. The zone-LP approach is a layout construction technique for the unequal-area departments and consists of two phases. In the first phase, a zoning algorithm is implemented to determine the relative positions between the departments. In this algorithm, for the sake of problem simplification and computational efficiency, each department is treated as a rectangle with an allowable aspect ratio and the area of the facility is assumed to be unbounded. In the second phase, by using the relative positions obtained in the first phase as input, a linear programming (LP) model is developed to identify the exact locations and dimensions of departments within the facility with specified sizes while satisfying their maximum aspect ratio requirement and the shape constraints. We also design a simulated annealing algorithm to improve the placing sequence. Finally, our computational results suggest that our proposed algorithm is efficient compared with the best existing approach in the literature. © 2016 PEI, University of Maribor. All rights reserved. Keywords: Facility layout problem Unequal area Zone-LP approach Simulated annealing *Corresponding author: zhanglm@nju.edu.cn (Zhang, L.M.) Article history: Received 10 August 2016 Revised 9 October 2016 Accepted 15 October 2016 1. Introduction As the competition among the global marketplaces increases rapidly, the provision of high-quality products or services at low cost to customers becomes more and more crucial for companies to capture the mass market. Facilities planning "determines how an activity's tangible fixed assets best support achieving the activity objectives" [1] and thus plays an important role in the cost reduction and productivity gain in industrial firms. Facility planning consists of facilities location, facility system design, facility layout design and handling system design. As a crucial part of facilities planning, facility layout design is to arrange departments interacting with each other within a facility so as to minimize the total material handling cost The layout design of a facility significantly impacts manufacturing cost and the productivity and efficiency of the system. According to [1], material handling cost ac- 259 Xiao, Zheng, Zhang, Kuo counts for 20 % to 50 % of the total manufacturing cost; an effective planning of facilities can reduce such cost by at least 10 % to 30 %. It is essential to design an efficient layout before manufacturing system design since future rearrangement of departments can result in a considerable expense. For these reasons, facility layout problem (FLP) has been studied extensively for over three decades [2-10]. There have been numerous studies adopting mathematical programming approaches to obtain optimal solutions for FLP. Montreuil [11] proposed the first mixed integer programming (MIP) model for solving the UA-FLP. However, their proposed MIP model can be solved for the problems with no more than six departments. Meller et al. [12] improved the formulation of the MIP model of Montreuil [11] by introducing a tightened department area constraint and several classes of valid inequalities. By their enhancement, the number of departments of the solvable problems increased to eight. Motivated by Meller et al. [12], Sherali et al. [13] further improved the accuracy for approximating the department areas by using a polyhedral outer approximation scheme. Moreover, Sherali et al. [13] investigated the effectiveness of the classes of valid inequalities proposed by Meller et al. [12] and reported that the MIP model was solved most efficiently by incorporating only two inequalities among them into the formulation. With this finding, the problems with up to nine departments could be optimally solved. Recently, Meller et al. [14] developed a new mathematical formulation for UA-FLP based on a sequence-pair representation. In this formulation, the feasible region of the problem was tightened such that the number of departments of the solvable problems has become eleven. Exact solution methods such as mathematical programming approaches are not applicable to obtain optimal solutions for large-scale problems since FLP is NP-hard [2]. Hence, heuristic approaches based on MIP models have been designed to construct layout solutions of high quality efficiently [15-17]. One of the well-known algorithms for solving FLP is based on the use of a slicing tree [18, 19], which is a binary tree used to form a subdivision of a rectangle by recursive computations. In the slicing tree, each inner node contains a guillotine cut operator (which cuts the area horizontally or vertically), and each leaf node represents a department. A layout can be generated according to a slicing tree by applying a set of guillotine cuts. Liu and Meller [20] developed a sequence-pair approach, where a pair of sequences is used to determine the relative locations between departments. By employing those relative locations, a reduced MIP model is utilized to construct the layout solutions. Both of above method are primarily proposed in the VLSI (very large Scale Integration) design. Similarly, Bozer and Wang [21] presented a graph-pair technique to fix the relative locations between departments and then a linear program (LP) is solved to obtain the layout solutions. This paper proposes a novel approach, which we call a combined zone-LP and simulated annealing (SA) algorithm, to solve the problem efficiently. The zone-LP approach is a layout construction technique for the UA-FLP and consists of two phases. In the first phase, the zone algorithm is implemented to determine the relative positions between departments. In this algorithm, for the sake of problem simplification and computational efficiency, each department is considered to be a rectangle with an allowable aspect ratio and the area of the facility is assumed to be unbounded. In the second phase, by using the relative positions obtained in the first phase as input, an LP is solved to determine the exact locations and dimensions of departments within the facility with specified sizes while satisfying their maximum aspect ratio requirement and the shape constraints. Furthermore, we also develop an SA algorithm to improve the solution quality. Finally, we conduct a computational study to examine the performance of our proposed algorithm and compare it with other existing solution approaches in the literature. Our paper is organized as follows. In the next section, we will describe the problem of UA-FLP. Section 3 will present our proposed solution methodology. In Section 4, we conduct computational experiments to study the performance of our proposed algorithm. Section 5 concludes our paper. 2. Problem description In UA-FLP, the areas of the departments are given and their dimensions are limited by a maximum aspect ratio (MAR), which is defined as the largest allowable ratio of the longest side to the 260 Advances in Production Engineering & Management 11(4) 2016 A combined zone-LP and simulated annealing algorithm for unequal-area facility layout problem shortest side of the department. The aspect ratio of a department ranges from 1/MAR to MAR. Furthermore, material flow quantities and the location of input/output (I/O) points are given. In the problem, the I/O points are assumed to be located at the centroids of departments. During the manufacturing process, materials are transferred between I/O points for each pair of departments. The objective is to minimize the total travel distance (TTD), which is the sum product of the distances between departments and their material flow quantities. For simplicity, the distance between departments is defined as the rectilinear distance between the input point and the output point. Table 1 provides an example of a setting of six departments with their required areas (am) and maximum aspect ratio (MARm) and the material flow quantities between each pair of the departments. Table 1 An example of a setting of six departments with their required areas (am) and maximum aspect ratio (MARm) and the material flow quantities between each pair of the departments m Material flow quantities 6 am MAR, 1 16 4 4 16 4 6 25 4 4 36 4 3 36 4 — 9 4 1 — 4 2 5 7 2 — — 4 3 3 3 — — — 2 2 4 7 6 — — — — — 3. Solution methodology It is well-known that UA-FLP is NP-hard. For this reason, we develop a heuristic algorithm to solve UA-FLP of large size. One complicating factor of UA-FLP is to determine the relative positions between the departments subject to the condition that their areas do not overlap. To circumvent this difficulty, we propose a solution methodology that integrates the zone-LP approach and an SA algorithm to obtain good-quality solutions. The zone-LP approach is a newly proposed layout construction technique based on an improved zone algorithm and the solution for an LP. When an initial solution is obtained from the previous stage, SA is employed to search within the solution space and acts to avoid the search procedure being trapped at a local optimum. 3.1 Zone-LP approach The zone-LP approach consists of two phases. In the first phase, for the simplification of the problem and more efficient computations, we consider the shape of each department as a rectangle with an allowable aspect ratio. As an example, as shown in Fig. 1, the six departments in Table 1 are treated as squares. Given the allowable aspect ratios of the departments, we develop an improved zone algorithm to determine their relative positions. The zone concept for UA-FLP is first introduced by Xiao et al. [22] to reduce its computational complexity and has been shown to be computationally efficient In the second phase, by using the relative positions obtained in the first phase as input, we solve an LP to obtain the exact locations of the departments within the facility, and their dimensions subject to their shape constraints are satisfied. 12 3 Fig. 1 Six departments to be located within the facility where each of them is of a fixed aspect ratio Advances in Production Engineering & Management 11(4) 2016 261 Xiao, Zheng, Zhang, Kuo 3.2 Improved zone algorithm The improved zone algorithm locates departments one by one according to a sequence such that a unique layout is generated. Given an initial sequence for locating the departments within the facility, say [1-4-2-5-3-6], we begin the algorithm by placing the first department in the sequence (i.e. Department 1) at the center of the facility. Then, four zones that are respectively on the left, above, on the right, and below this department can be identified. As an example, the four zones are shaded in the graphs of Fig. 2. The next step is to determine the zone that the following department in the sequence (i.e. Department 4) to be located. We first determine the optimal locations of this department within the four zones by using a median point (MP) method. The department will then be located at the best zone that results in the minimum total distance. Once the department is located within the facility, the set of zones will be updated for the determination of the location of the next department in the sequence. That procedure is repeated until all departments are located. For more detail about zone updating procedure, we refer the reader to Xiao etal. [22]. (a) (b) Fig 2 The four zones (shaded in grey) for the possible location of the department in our example To generalize the expressions of the zones, for i = 1, ...,N, let Ft be the i-th department in the sequence and let Z(i) = {z1l,.,zjl,.,zj1} be the set of zones that Ft can be located, where j is the index of zone. Note that Z^ consists of only one zone that is the entire rectangle facility. For each zone of Fj, Zjl is represented by where(xlj1,ylj1) and (xlj2,ylj2) are the coordinates of the bottom left corner and the top right corner of the zone, respectively. For step of determining the location of the department, the MP method is applied to obtain the optimal position of Ft within each zone in Z^. In this MP method, an ideal location of the cen-troid of Fi is first identified in such a way that the sum of rectilinear distances weighted by the flow quantities between Ft and the departments that have previously been located, denoted by the set B, is minimized. The objective function is formulated as follows: Minimize fmn(\xi ~xn\ + \yj _yJ) (1) neB To determine the optimal location, (x*,y*), the objective function Eq. 1 is reformulated as Equations Eq. 2 and Eq. 3. Therefore, the values of x* and y* can be calculated separately. In Equations Eq. 2 and Eq. 3, TTDi(x) and TTDi(y) are piece wise linear and convex functions, where all xn (yn), for neB, are arranged in increasing order, i.e. x±'< Xt'< — W/2. Hence, the MP is x3' = 0 and x* = x3' = 0. Similarly, TTDs(y) = 2|y5-0| + 2|y5-0| + 2|y5-0| + 4|ys-4|. And Wy(1) = 2, Wy(2) = 2 + 2 = 4, Wy(3) = 2 + 2+ 2 = 6> W/2 . Hence, the MP in y -direction is y3' = 0 and yl = y-i' = 0. Therefore, the optimal location of Department 3 is (0, 0), which is marked by the dashed rectangle in Fig. 3. |B| ,|B| Fig. 3 Optimal location of Department 3 (represented by a dotted square) Nevertheless, this optimal location determined for Department 3, in fact, is infeasible because it overlaps with the other departments as shown in in Fig. 3. In this case, we have to determine a new location for Department 3 to continue our UA-FLP procedure. To this end, we determine the location in each zone for Department 3 such that the location is closest to the previously determined optimal location in both x-direction and y-direction. The rationale is that such position can minimize TTD for each zone. For example, the best positions for Department 3 in those considered zones are shown in Fig. 4(a)-4(f), which are respectively the closest locations to the optimal location that was previously determined. The department is to be located at the best position that attains the minimum TTDij among those zones. In this example, Department 3 has the same TTD^ at the positions as shown in Fig. 4(b) and Fig. 4(d). Thus, this tie is broken by randomly picking one of the two positions for locating Department 3. The location of the last department in the sequence of this example is also determined in the same way. Fig. 5 shows the final layout of the facility. Advances in Production Engineering & Management 11(4) 2016 263 Xiao, Zheng, Zhang, Kuo (a) y (b) 4 2 4 5 1 -5 0 7 (e) (d) (e) (f) Fig. 4 Closest location in each zone to the optimal location y -5 Fig. 5 Layout of the six departments of fixed shapes To determine if the department can be located at the intended position, we also have to consider the dimension of each zone zf which can be measured by wjl = xlj2 and h} =ylj2 ~y)i, where wjl and h} are respectively the width and the height of the zone. Let dwm and dhm be the width and the height of the department to be located, respectively. If a zone satisfies the condi- tions that Wjl >dwm and h- >dhm, it is feasible to locate Department m in this zone; otherwise, check if the subsequent zone can accommodate Department m. The procedure of the modified zone algorithm is summarized as follows: 264 Advances in Production Engineering & Management 11(4) 2016 A combined zone-LP and simulated annealing algorithm for unequal-area facility layout problem Step 1. Set i = 1. The centroid of F( is placed at (0, 0) in the coordinate system. Step 2. Set i = i + 1. Create the zone set, Z^ = [zjl\j = 1,...,/}, by using the zone updating procedure. Determine the optimal location, (x*, y*), for F(. Let TTD* be the minimum TTDij among the zones contained in the zone set. Set j = 1, TTD* = +<». Step 3. If xlJi + dwm/2 < x* < xlJ2-dwm/2 and ylJi + dhm/2 < y* < y)2 -dhm/2 (i.e. F^ with its centroid located at the optimal location can fit entirely within zone zlj), go to Step 3.1; otherwise, go to Step 3.2. Step 3.1 F( is optimally located at (x*, y*). Go to Step 5. Step 3.2 If wlj ^dwm and wlj ^ dhm, locate Di at the closest possible position in Zjl to the optimal location. Calculate TTDtj. If TTDtj 2 C penalty for the violation of the floor boundary condition, c= XZfmn(w+H) m n>m wm,hm Half of width and height of Department m vxm, vym Amount of violation of the floor boundary condition for Department m in the x-direction and the y-direction , respectively Our LP model to determine the dimensions of the departments is derived from the MIP formulated by Sherali et al. [13] and presented below: Minimize TTD = ^ ^ fmni^-mn + + Vm + ^ (4) m n>m m m subject to Advances in Production Engineering & Management 11(4) 2016 265 Xiao, Zheng, Zhang, Kuo ^mn —^ji ^Tn ^ ^ (6) rfmn ^yin-Vil V™ < n (7) dmn ^Yn-ym Vm0 Vm (17) ^mrn^inn — 0 Vm < n (18) Objective function Eq. 4 consists of two terms: TTD between departments and the penalty for violation of the floor boundary condition. Constraints Eq. 5 to Eq. 8 determine the TTD between departments. Constraints Eq. 9 and Eq. 10 measure the violation of the floor boundary condition (in x-direction and y-direction, respectively), if it is infeasible to locate all the departments within the facility. Constraints Eq. 11 and Eq. 12 impose the bounds on the width and height for each department, where the upper bounds ub^ = min{W, ^Jamam}, ub^ = min{7/, ^Jamam}, and lower bounds lb^ = am/ublb^ = am/ubxm. Constraints Eq. 13 and Eq. 14 ensure that the departments do not overlap with each other. Constraints Eq. 15 is the polyhedral outer approximation for the area function, at = 4wmhm, with the number of tangential support points equal to A. The value of the tangential support points x is given in Eq. (16). As an example, with the relative positions of the departments obtained by the improved zone algorithm (as shown in Fig. 5), the final facility layout for locating the six departments with the exact dimensions can now be obtained by the LP. Their positions and dimensions are shown in Fig. 6. y 266 Fig. 6 Layout solution of the six departments with their dimensions determined Advances in Production Engineering & Management 11(4) 2016 A combined zone-LP and simulated annealing algorithm for unequal-area facility layout problem 3.4 SA algorithm With the facility layout solution obtained by the zone-LP approach, we implement an SA algorithm to improve the quality of the solution. The SA algorithm that we adopt is the same as the procedure proposed in Xiao et al. [22]. For the detail of the algorithm, we refer the reader to Xiao et al. [22]. 4. Computational experiments In this section, we conduct computational experiments to evaluate the performance of our proposed solution methodology for UA-FLP - an integrated method of zone-LP approach and SA algorithm. We use a set of widely tested instances in the literature to examine the efficiency of the solution methodology. There are five instances of different problem sizes included in the computational experiments: Problems 07 (with 7 departments), 08 (with 8 departments) and 09 (with 9 departments) from Meller et al. [12] and Two largest instances SC30 (with 30 departments) and SC35 (with 35 departments) from Liu and Meller [20]. The combined zone-LP and SA algorithm is run ten times for each computational instance. In the first phase of zone-LP algorithm, the dimensions of departments are fixed with the MAR because all departments tend to have a narrower rectangular shape to get a lower TTD. The parameter combination used in SA was T = 200 R = 0.8, and L = 1000 from Xiao et al. [21]. The results on the computational performance of our proposed approach are shown in Table 2. We present the averages and the standard deviations of the costs to examine the robustness of our proposed method. The computational results suggest that the proposed algorithm appears to be quite robust. For the three small-sized problems (07, 08, 09), the standard deviation is quite small, ranging from 0 to 0.57 (less than 0.5 % of the average). Although the magnitude of the standard deviation becomes large for the large-sized problems (SC30, SC35), its relative value is still small (less than 3 % of the average). For the computational efficiency, our proposed solution methodology can obtain the final solution within a reasonable time frame (less than 1 minute for the small-sized problems and less than 4 minutes for the large-sized problems). Given that facility layout planning is not a daily practice (is usually determined for several years), this computational time is negligible, and our method is suitable for the application in practice. Table 2 Qbjective values and the computing time obtained by the proposed algorithm Problem Best Mean Worst Standard deviation Average time (s) O7 89.25 89.25 89.25 0 33.13 O8 185.00 185.30 186.00 0.46 36.93 O9 185.00 185.45 186.5 0.57 41.80 SC30 3441.57 3663.21 3792.71 103.50 180.27 SC35 3347.94 3423.70 3555.89 69.06 225.94 Table 3 0bjective values, dimensions of the facilities, space utilizations resulting from the best solutions found by GRAPH [21] and our proposed algorithm GRAPH[21] The proposed algorithm Diff. in TTD (%) Problem TTD Layout dimension (WxH) Space utilization (%) TTD Layout dimension Space utilization (WxH) (%) O7 115.93 96.00 89.25 12x13.5 68.52 23.01 O8 239.00 96.00 185.00 12x16.5 74.24 22.59 O9 227.10 96.00 185.00 15x16.5 63.03 18.54 SC30 3601.20 15x12 90.56 3441.57 12.75x21.59 59.21 4.43 SC35 3351.12 16x15 80.00 3347.94 15.14x24.7 51.34 0.09 We also compare our proposed solution methodology with other existing algorithms in the literature and use them as benchmarks. We first compare the solution qualities obtained by our approach and the GRAPH algorithm [21], which attains the best known solutions for the computational instances. As shown in Table 3, our proposed outperforms GRAPH [21] in terms of TTD for all computational instances, but the space utilizations resulting from our algorithm are less than those from GRAPH. In our UA-FLPs, departments are arranged within an open-field floor Advances in Production Engineering & Management 11(4) 2016 267 Xiao, Zheng, Zhang, Kuo (as shown in Fig. 4 and Fig. 5) and the loss of consideration on space utilization in the objective function of the proposed model mainly causes the poor performance on space utilization. Therefore, a multi-objective optimization considering minimization of TTD and maximization of space utilization may be our future research topic. For reference, the best layouts solutions obtained by our proposed approach are shown in Appendix 1. The proposed algorithm also has an advantage of a shorter computing time due to its capability of searching in the restricted solution spaced. We conduct another study to examine the efficiency of our proposed methodology. In addition to GRAPH [21], we also include SEQUENCE [20] for the comparison of computational efficiencies of the different approaches. Computational results in Table 4 suggest that our algorithm takes a significantly shorter time to solve the UA-FLP. More importantly, the computing time appears to grow almost linearly as the number of departments increases. Compared with the approaches that appear to have an exponential computing time in the number of departments, our approach is apparently more suitable for problems of a larger size. The linear increase of computing time may be attributed to the proposed zone-LP algorithm which can construct the layout solution effectively. In each iteration of SA, the number of zones produced for putting departments by MP method is O (N) according to the zone algorithm. The proposed LP model is just used to determine the dimensions of departments. In general, the computing time is approximately proportional to the number of generated zones in each SA iteration and thus grows almost linearly as the number of departments increases. This further demonstrates the value of our proposed approach in advancing the computational performance for solving UA-FLP. Table 4 Average computing times of the SEQUENCE [20], GRAPH [21], and our proposed algorithm Instance Computing time (s) SEQUENCE [20] GRAPH [21] Our proposed algorithm O7 1644.0 228.0 33.13 O8 3056.0 390.0 36.93 O9 3879.0 222.0 41.80 SC30 7282.8 2442.0 180.27 SC35 9590.4 4728.0 225.94 5. Conclusion and future research This paper deals with UA-FLP and proposes a novel approach, which we call a combined zone-LP and SA algorithm, for solving large-sized UA-FLP. The zone-LP algorithm is a new layout construction method and consists of two phases. In the first phase, the shapes of departments are fixed to a rectangle with an allowable aspect ratio. Those departments of fixed shapes are then placed sequentially by using an improved zone algorithm, where an MP method is proposed to locating departments. By using relative positions of the departments obtained in this first phase, an LP is formulated to determine the exact locations and dimensions of departments. We also implement an SA algorithm to search the sequences of locating the departments within the facility and to prevent the search procedure from getting trapped at a local optimum. Computational experiments indicate that our proposed combined zone-LP and SA has a reasonably good performance, compared with existing algorithms in the literature on widely tested computational instances. We also note that there can be future directions extended from this research. In reality, departments of irregular shapes are commonplace; a future study may consider the facility layout problem with the consideration of departments with irregular shapes. Moreover, the results of computing experiments suggest that our proposed solution methodology may generate solutions leading to a low space utilization. For this case, we may adopt a multi-objective optimization approach whose goals are to minimize the material handling cost and to maximize space utilization. The trade-off between these two performance measures would also be an interesting research topic for investigation in the future. 268 Advances in Production Engineering & Management 11(4) 2016 A combined zone-LP and simulated annealing algorithm for unequal-area facility layout problem Acknowledgement The research has been supported by the National Natural Science Foundation of China (Grant No. 71501090, 71501093), the Natural Science Foundation of Jiangsu Higher Education Institutions of China (Grant No. 14KJD410001), the Natural Science Foundation of Jiangsu Province (Grant No. BK20150566), and the General Research Fund of Research Grants Council of Hong Kong (Grant No. 414313). References [1] Tompkins, J.A., White, J.A., Bozer, Y.A., Tanchoco, J.M.A. (2010). Facilities Planning, 4th edition, John Wiley & Sons, New York, USA. [2] Kusiak, A., Heragu, S.S. (1987). The facility layout problem, European Journal of Operational Research, Vol. 29, No. 3, 229-251, doi: 10.1016/0377-2217(87)90238-4. [3] Meller, R.D., Gau, K.Y. (1996). The facility layout problem: Recent and emerging trends and perspectives, Journal of Manufacturing Systems, Vol. 15, No. 5, 351-366, doi: 10.1016/0278-6125(96)84198-7. [4] Singh, S.P., Sharma, R.R.K. (2006). A review of different approaches to the facility layout problems, The International Journal of Advanced Manufacturing Technology, Vol. 30, No. 5-6,425-433, doi: 10.1007/s00170-005-0087-9. [5] Drira, A., Pierreval, H., Hajri-Gabouj, S. (2007). Facility layout problems: A survey, Annual Review in Control, Vol. 31, No. 2, 255-267, doi: 10.1016/j.arcontrol.2007.04.001. [6] Kulturel-Konak, S., Konak, A. (2011). Unequal area flexible bay facility layout using ant colony optimisation, International Journal of Production Research, Vol. 49, No.7, 1877-1902, doi: 10.1080/00207541003614371. [7] Jolai, F., Tavakkoli-Moghaddam, R., Taghipour, M. (2012). A multi-objective particle swarm optimisation algorithm for unequal sized dynamic facility layout problem with pickup/drop-off locations, International Journal of Production Research, Vol. 50, No. 15, 4279-4293, doi: 10.1080/00207543.2011.613863. [8] Navidi, H., Bashiri, M., Bidgoli, M.M. (2012). A heuristic approach on the facility layout problem based on game theory, International Journal of Production Research, Vol. 50, No. 6, 1512-1527, doi: 10.1080/00207543.2010. 550638. [9] Ficko, M., Palcic, I. (2013). Designing a layout using the modified triangle method, and genetic algorithms, International Journal of Simulation Modelling, Vol. 12, No. 4, 237-251, doi: 10.2507/IISIMM12(4)3.244. [10] Ficko, M., Brezovnik, S., Klancnik, S., Balic, J., Brezocnik, M., Pahole, I. (2010). Intelligent design of an unconstrained layout for a flexible manufacturing system, Neurocomputing, Vol. 73, No. 4-6, 639-647, doi: 10.1016/ j.neucom.2009.06.019. [11] Montreuil, B., 1991. A modelling framework for integrating layout design and flow network design, Progress in Materials Handling Research, Vol. 2, 95-115, doi: 10.1007/978-3-642-84356-3 8. [12] Meller, R.D., Narayanan, V., Vance, P.H. (1998). Optimal facility layout design, Operations Research Letters, Vol. 23, No. 3-5, 117-127, doi: 10.1016/S0167-6377(98)00024-8. [13] Sherail, H.D., Fraticelli, B.M.P., Meller, R.D. (2003). Enhanced model formulations for optimal facility layout, Operations Research, Vol. 51, No. 4, 629-644, doi: 10.1287/opre.51.4.629.16096. [14] Meller, R.D., Chen, W., Sherali, H.D. (2007). Applying the sequence-pair representation to optimal facility layout design, Operations Research Letters, Vol. 35, No. 5, 651-659, doi: 10.1016/j.orl.2006.10.007. [15] Kim, J.G., Kim, Y.D. (2000). Layout planning for facilities with fixed shapes and input and output points, International Journal of Production Research, Vol. 38, No. 18, 4635-4653, doi: 10.1080/00207540050205550. [16] Das, S.K. (1993). A facility layout method for flexible manufacturing systems, International Journal of Production Research, Vol. 31, No. 2, 279-297, doi: 10.1080/00207549308956725. [17] Rajasekharan, M., Peters, B.A., Yang, T. (1998). A genetic algorithm for facility layout design in flexible manufacturing systems, International Journal of Production Research, Vol. 36, No. 1, 95-110, doi: 10.1080/00207549 8193958. [18] Scholz, D., Petrick, A., Domschke, W. (2009). STaTS: A Slicing Tree and Tabu Search based heuristic for unequal area faciltiy layout problem, European Journal of Operational Research, Vol. 197, No. 1, 166-178, doi: 10.1016/ j.ejor.2008.06.028. [19] Komarudin, Wong, K.Y. (2010). Applying ant system for solving unequal area facility layout problems, European Journal of Operational Research, Vol. 202, No. 3, 730-746, doi: 10.1016/j.ejor.2009.06.016. [20] Liu, Q., Meller, R.D. (2007). A sequence-pair representation and MIP model based heuristic for the facility layout problem with rectangular departments, IIE Transactions, Vol. 39, No. 4, 337-394, doi: 10.1080/07408170600 844108. [21] Bozer, Y.A., Wang, C.T. (2012). A graph-pair representation and MIP-model-based heuristic for the unequal-area facility layout problem, European Journal of Operational Research, Vol. 218, No. 2, 382-391, doi: 10.1016/ j.ejor.2011.10.052. [22] Xiao, Y., Seo, Y., Seo, M. (2013). A two-step heuristic algorithm for layout design of unequal-sized facilities with input/output points, International Journal of Production Research, Vol. 51, No. 14, 4200-4222, doi: 10.1080/ 00207543.2012.752589. Advances in Production Engineering & Management 11(4) 2016 269 Xiao, Zheng, Zhang, Kuo Appendix 1 The best layouts obtained by our proposed algorithm for the experimental instances from Meller et al. [10] and Liu and Meller [19]. Layout with seven departments Layout with eight departments Layout with nine departments 25 r 15 a 20 19 23 zn ii i 13 3 29 27 26 0 2 4 6 8 10 12 14 Layout with 30 departments 16 1D "26 2« 19 £IE 15 13 14 32 12 H 20 30 29 22 34 27 35 ! 4 6 8 10 12 14 Layout with 35 departments 270 Advances in Production Engineering & Management 11(4) 2016 APEM jowatal Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 | pp 271-286 http://dx.doi.Org/10.14743/apem2016.4.226 ISSN 1854-6250 Journal home: apem-journal.org Original scientific paper A new multi-objective Jaya algorithm for optimization of modern machining processes Rao, R.V.a*, Rai, D.P.a, Ramkumar, J.b, Balic, J.c aDepartment of Mechanical Engineering, Sardar Vallabhbhai National Institute of Technology, Surat, India bDepartment of Mechanical Engineering, Indian Institute of Technology, Kanpur, India Production Engineering Institute, Faculty of Mechanical Engineering, University of Maribor, Slovenia A B S T R A C T A R T I C L E I N F O In this work, the multi-objective optimization aspects of plasma arc machining (PAM), electro-discharge machining (EDM), and micro electro-discharge machining (|-EDM) processes are considered. Experiments are performed and actual experimental data is used to develop regression models for the considered machining processes. A posteriori version of Jaya algorithm (MO-Jaya algorithm) is proposed to solve the multi-objective optimization models in a single simulation run. The PAM, EDM and |-EDM processes are optimized using MO-Jaya algorithm and a set of Pareto-efficient solutions is obtained for each of the considered machining processes and the same is reported in this work. This Pareto optimal set of solutions will provide flexibility to the process planner to choose the best setting of parameters depending on the application. The aim of this work is to demonstrate the performance of MO-Jaya algorithm and to show its effectiveness in solving the multi-objective optimization problems of machining processes. © 2016 PEI, University of Maribor. All rights reserved. Keywords: Plasma arc machining Electro-discharge machining Micro-electro-discharge machining Multi-objective optimization Jaya algorithm Posteriori approach Sustainability *Corresponding author: ravipudirao@gmail.com (Rao, R.V.) Article history: Received 10 November 2016 Revised 25 November 2016 Accepted 1 December 2016 1. Introduction In order to survive in a fierce market scenario manufacturing industries are required to maintain high quality standards, produce at lowest cost, increase production rate, conserve resources and at the same time minimize the environmental impact of the processes they use. Machine tools are major pillars of any manufacturing system and are used on a large scale for processing of materials. However, machining processes are characterized by high energy consumption, high tool wear rate, poor surface quality and generation of large scale waste products in the form of used lubricants, coolants, dielectric or electrolytic fluids, chips and debris of tool or workpiece materials, etc. Thus, for success of any manufacturing system in terms of economy and to reduce its impact on the ecology it is crucial to improve the efficiency of these machine tools. Furthermore, in order to improve the sustainability of the process it is imminent that the machines are operated as efficiently as possible. The performance of any machining process extensively depends upon the choice of process parameters. Therefore, for best performance from any machining process it is important to set the process parameters optimally. In order to determine the optimal setting of process parameters it is important to map the relationship between input and output parameters. De Wolf et al. [1] investigated the effect of process parameters on material removal rate, electrode wear rate and surface finish in EDM process. Aich and Banerjee [2] applied teaching learning based opti- 271 Rao, Rai, Ramkumar, Balic mization procedure for the development of support vector machine learned EDM process and its pseudo Pareto optimization. Zhang et al. [3] enumerated and characterized 128 scenarios in sustainable machining operation involving 7 objectives including energy, cost, time, power, cutting force, tool life and surface finish. Gupta et al. [4] presented the results of optimization of machining parameters and cutting fluids during nano-fluid based minimum quantity lubrication turning of titanium alloy by using particle swarm optimization and bacteria foraging optimization techniques. Researchers have also applied a number of numerical and metaheuristic optimization algorithms for optimal setting of machining process parameters [5-13]. The metaheuristic optimization algorithms are mostly inspired by the theory of evolution or of behavior of a swarm. All evolutionary algorithms or swarm based algorithms require tuning of parameters like population size, number of iterations, elite size, etc. In addition, different algorithms require their own algorithm-specific parameters. The improper tuning of algorithm-specific parameters adversely affects the performance of these algorithms. In addition, the tuning of population size and number of iterations is also required. Rao [14] proposed the Jaya algorithm which algorithm-specific parameter-less algorithm. The performance of Jaya algorithm has already been tested on a number of unconstrained and constrained benchmark functions and engineering optimization problems. For more details about the algorithm, the readers may refer to https://sites.google.com/site/jayaalgorithm. The Jaya algorithm is simple in implementation as a solution is updated only in a single phase using a single equation. However, the multi-objective version of Jaya algorithm is not yet developed. In the case of machining processes due to co-existence of multiple performance criteria there is a need to formulate and solve multi-objective optimization problems (MOOP). A priori approach such as normalized weighted sum approach, epsilon constraint method, etc. require assigning the weights of importance to the objectives before simulation run of the algorithm. Further, it is required to run the algorithm independently for each set of weights to obtain distinct solutions. A posteriori approach does not require assigning weights of importance to the objectives in advance. This approach provides a set of Pareto-efficient solutions for a MOOP in a single run of simulation. The process planner can then select one out of the set of Pareto-efficient solutions based on the order of importance of objectives. Thus, in this work a parameter-less posteriori multi-objective version of Jaya algorithm is named as multi-objective Jaya (MO-Jaya) algorithm is proposed and the MOOPs of three modern machining processes namely plasma arc machining (PAM), electro-discharge machining (EDM), and micro electro-discharge machining (^-EDM) are solved using MO-Jaya algorithm. The Jaya and MO-Jaya algorithms are described in following sections. 2. The Jaya algorithm In the Jaya algorithm P initial solutions are randomly generated obeying the upper and lower bounds of the process variables. Thereafter, each variable of every solution is stochastically updated using Eq. 1. The best solution is the one with maximum fitness (i.e. best value of objective function) and the worst solution is the one with lowest fitness (i.e. worst value of objective function). Op+i,q,r Op,q,r 1 p,q,best abs(^p,q,r)) ^p,q,2 p,q,worst abs(^p,q,r)) (1) Here best and worst represent the index of the best and worst solutions among the population. p, q, r are the index of iteration, variable, and candidate solution. Op, q, r means the q-th variable of r-th candidate solution in p-th iteration. ap,q,t and ap,q,2 are numbers generated randomly in the range of [0, 1]. The random numbers ap,q,1 and ap,q,2 act as scaling factors and ensure exploration. The absolute value of the variable (instead of a signed value) also ensures exploration. Fig. 1 gives the flowchart for Jaya algorithm. 272 Advances in Production Engineering & Management 11(4) 2016 A new multi-objective Jaya algorithm for optimization of modern machining processes Fig. 1 Flowchart of Jaya algorithm 3. The multi-objective Jaya algorithm The MO-Jaya algorithm is a posteriori version of Jaya algorithm for solving MOOPs. The solutions in the MO-Jaya algorithm are updated in the similar manner as in the Jaya algorithm based on Eq. 1. In the interest of handling problems in which more than one objective co-exist the MO-Jaya algorithm is embedded with dominance ranking approach and crowding distance evaluation approach. The MO-Jaya algorithm is a posteriori version of Jaya algorithm for solving MOOPs. The solutions in the MO-Jaya algorithm are updated in the similar manner as in the Jaya algorithm based on Eq. 1. In the interest of handling problems in which more than one objective co-exist the MO-Jaya algorithm is embedded with dominance ranking approach and crowding distance evaluation approach [12]. In the MO-Jaya algorithm, the superiority among the solutions is decided according to the non-dominance rank and value of the density estimation parameter i.e. crowding distance (%). The solution with highest rank (rank = 1) and largest value of % is chosen as the best solution. On the other hand the solution with the lowest rank and lowest value of % is selected as the worst solution. Such a selection scheme is adopted so that solution in less populous region of the objective space may guide the search process. Once the best and worst solutions are selected, the solutions are updated based on the Eq. 1. After all the solutions are updated, the updated solutions are combined with the initial population to so that a set of 2P solutions (where P is the size of initial population) is formed. These solutions are again ranked and the % value for every solution is computed. Based on the new ranking and % value P good solutions are chosen. The flowchart of MO-Jaya algorithm is given in Fig. 2. For every candidate solution the MO-Jaya algorithm evaluates the objective function only once in each iteration. Therefore, the total no. of function evaluations required by MO-Jaya algorithm = population size x no. of iterations. However, when the algorithm is run more than once, then the number of function evaluations is to be calculated as: no. of function evaluations = no. of runs x population size x number of iterations. The methodology used for ranking of solutions, computing the crowding distance and crowding comparison operator are described in the following sub-sections. Advances in Production Engineering & Management 11(4) 2016 273 Rao, Rai, Ramkumar, Balic Fig. 2 Flowchart of MO-Jaya algorithm 3.1 Ranking methodology The approach used for ranking of solutions is based on the non-dominance relation between solutions and is described as follows. In an M objective optimization problem, P is the set of solutions to be sorted and n = | P|. Domination: A solution xi is said to dominate another solution x2 if and only iff (xi) < f (x2) for all 1 < i < M andf (x1) Model - 4 Y?" Model - 5 max | I | min Fig. 6 Total displacement distribution of the steel bearing components Advances in Production Engineering & Management 11(4) 2016 293 Özkal, Cakir, Arkun Modal Snow Wind o -a o S CS I "o -a o S CO I -a o S i o -a o S i -a o S max min Fig. 7 Total displacement distribution of the all structural components 6. Determination of the optimum design A structural designer should consider the effective use of materials along with the safety limit and esthetical semblance objectives. Although this point has been investigated in detail by several researchers over the last few decades, manufacturers generally put the semblance forward for the purpose of commercial concerns and push the structural performance of their designs into background. If a technique is sought in order to answer the question "Which one is the best?", a performance qualification method is needed to be constituted subsequent to the definition of all the design criterions [21]. In this study, five types of popular carport designs covering equal areas were evaluated. Because expected consumer benefit equality is provided for all of the models, there is no need to consider this factor. Furthermore, it should be noted that each product, including carport structures, must be optimized not only based on material concerns or structural behaviour, but also by considering the easy-for-manufacture and easy-for-assembly paradigms. The most suitable design is usually a compromise among the above mentioned requirements. Because current study deals with the rough design at the pre-production stage, it has also been neglected for the aim of obtaining optimum result. For instance, a stability analysis incorporating structural assembly details that was studied by Manifold [22] or a geotechnical investigation that was studied by Hrestak et al. [23] could be integrated as a post-verification stage to the design process in order to make sure every constraint is satisfied before the production. 294 Advances in Production Engineering & Management 11(4) 2016 Finite element method for optimum design selection of carport structures under multiple load cases Design objective for a carport structure could be summarized as satisfying material efficiency objective while maintaining a rigid behaviour against various possible load cases. Although uniform stress levels for the whole structure will assist material efficiency, size and shape dissimilarity of finite elements prevent an accurate calculation. Hence, a new modified performance index formulation, based generally on the nodal displacement levels, should be built according to the following problem definition: minimize W = £e=iP Ve subject to umax1 umax U >U* ^ ^ >1 avg — avg ,,* — J- uavg where W is the weight of the structure, p is the material density, Ve is the volume of the eth element, t is the total number of elements, umax is the absolute value of the maximum nodal displacement, uavg is the absolute value of the average nodal displacement value while u*^iax and u*avg are the upper and lower bound limits of the displacement constraints, respectively. For the determination of the optimality of carport designs, a performance index, which could be applied in a general scope and considers all the load cases, should be developed step by step. In order to constitute a rigid structure, structural performance level based on the maximum displacement value of the zth model for the jth load case is firstly defined. Xi 7 = 1 U max Secondly, the average nodal displacement value should be added to the formulation for the purpose of assisting material efficiency qualification. XÍ' ZI umax uavg \ i/J u* I (2) 7 = 1 \Umax av3 The mathematical definition of structural performance level is finalized by implementing the weight value as follows. 1 n ( * \ _ 1 \ | U^ax uavg \ 1 7 = 1 \ 9 max/ However, there is a need to define a reference performance level in order to compare the designs with each other. Considering a certain number of models are evaluated in this study, this reference level could be formulated as the mean performance level of all models. 1 m ¿=1 (4) Finally, a performance index formulation is generated by the proportion of the zth model's and reference performance level values. It is named as PIrd because of the aim to attain a rigid design while minimizing the material weight of the structure; in other words, while maximizing the material efficiency. PIrd — ~ — _ Xref 1 _Ly« i Uxijk (11) i = 0 j = 0 In Eq. 10, the waiting time and operating time are given weights, and in Eq. 11, the waiting time is given the weight. Min(Max(Tei)),ieV (12) Max(Tei) >TSi — Tfi + WSi/BVj (13) (TSi -Tf + WSi/BVj)-(Tsi+1 -Tfi+1 + Wsi+1/BVj+1)>0 (14) Eq. 12 minimizes the makespan. 4. MPPSO for BAP and computational experiments 4.1 Algorithm design MPPSO is developed to solve the BAP for a bulk cargo port. (1) Position coding Each vessel has information about the job berths and job sequence, so each particle has two-dimensional information about the job berths and job sequence. Job berths are affected by a vessel's own properties and berth conditions. Thus, we need to obtain the berths arrays firstly. For example, if berths 1, 3, and 6 are available for each vessel, the berths array is [1 3 6 0 0 0] (assuming six berths). The job sequence of particles can be explained by the continuous value; the smaller the value is, the earlier the vessel will obtain service (more details in Table 1). Table 1 Example of position coding Dimension 1 2 3 Berths [1 3 6 0 0 01 [1 3 0 0 0 01 [6 0 0 0 0 0 1 Sequence 2 0.6173 2 0.3214 1 0.4142 In Table 1, we can obtain the job berths and job sequence through the position decoding; that is, vessels 1 and 2 obtain service in berth 3, where berth 3 provides service for vessel 2 firstly, and vessel 3 obtains service in berth 6. (2) Position renewal The parameter setting of MPPSO is as follows: • T: maximum iterations • N: number of particle swarms • D: dimension • V;d(t): velocity of the i-th particle with the d-th dimension at time t • xid(t): position of the i-th particle with the d-th dimension at time t • gd(t): best positions discovered by all particles at time t or earlier, with d-th dimension • ph: number of phases • pcf: phase change frequency • g: number of groups within each phase • Cv , Cx, and Cg: coefficient value in each group within each phase • sl: sub-length of the dimension • VC: velocity change variable Therefore, the renewal process of velocity and position is expressed as follows: vidi(.t + 1) = Cvvidl(t) + Cxxidl(t) + Cggdl(t) (15) Xidi (t + 1) = *idi(0 + Vidi(t + 1) (16) Advances in Production Engineering & Management 11(4) 2016 303 Tang, Gong, Liu, Zhang vid2(t + 1) = cvvid2(t) + Cxxid2(t) + Cggd2(t) (17) xid2(t +1) = xid2(t)+ vid2(t +1) The pseudo-code for MPPSO is shown in Table 2: (18) Table 2 Pseudo-code for MPPSO Step 0 read data: Vn and Bn Step 1 get the berths arrays based on data Step 2 T D N ph pcf g VC Step 3 initialize the velocity and position for particles, get the best of positions Pgx and the best fitness value Pg Step 4 iteration for t = 1 : T if t is a multiple of VC reinitiate velocity determine its phase for i = 1 : N determine its group choose sl randomly from [1,min(10,n)] sl = roundn(rand*(min(10,D)-1)+1,0); deal with each dimension of particle d = 1; initialize the current dimension dimension limited while ( d <= D ) Cache the initial position of particle temp (1,: ,:) = x(i,: ,:); deal with the sub dimension of particle for j = 0 : sl updating dimension d_temp = d + j; if d_temp >D,break Step 4.1 updating the berth get the coefficient of velocity formula update the velocity of berth based on formula(12) update the position of berth based on formula(13) if out of berth array scope, reinitiate position Step 4.2 updating the sequence get the coefficient of velocity formula update the velocity of vessel based on formula(14) update the position of vessel based on formula(15) end judge the fitness if fitness(temp) < fitness(x(i,: ,:)) accept the updating x(i,: ,:) = temp(1,: ,:); end update the best of positions in particle swarm if fitness(x(i,: ,:)) max_finish_time max_finish_time = berths_ships(i,7); end 0 200 W0 600 800 1000 1200 U00 1600 1800 2000 geneiaion (a) 0 200 <00 «10 B00 10TO 1200 "00 1600 1BTO 2000 generation 200 ÍOO 600 BOO 1000 1200 M00 1600 1930 2030 0 200 4O0 600 600 1000 1200 1100 1600 1800 2000 genefaliwi M fd) Fig. 1 Iteration process and fitness value Table 10 MPPSO results Model Berth ID Vessel ID Starting time Waiting time Finishing time Operating time 1 2014-08-13 01:00 0.00 2014-08-13 15:47 0.62 1 4 2014-08-13 15:47 0.16 2014-08-14 04:40 0.54 12 2014-08-14 04:40 0.40 2014-08-14 17:00 0.51 13 2014-08-13 00:00 0.00 2014-08-13 15:00 0.63 2 10 2014-08-13 15:00 0.42 2014-08-14 06:51 0.66 9 2014-08-14 06:51 0.33 2014-08-14 20:24 0.56 3 3 2014-08-13 09:30 0.00 2014-08-13 10:15 0.03 4 5 2014-08-13 05:00 0.00 2014-08-14 00:45 0.82 Advances in Production Engineering & Management 11(4) 2016 307 1 Tang, Gong, Liu, Zhang Table 10 MPPSO results (continuation) 7 2014-08-13 07:00 0.00 2014-08-13 16:56 0.41 11 2014-08-13 17:00 0.00 2014-08-14 03:54 0.45 6 2014-08-14 03:54 0.33 2014-08-14 07:30 0.15 8 2014-08-14 07:30 0.56 2014-08-14 18:44 0.47 14 2014-08-13 00:00 0.00 2014-08-13 16:36 0.69 2 2014-08-13 19:00 0.00 2014-08-14 14:15 0.80 1 2014-08-13 01:00 0.00 2014-08-13 15:47 0.62 1 11 2014-08-13 17:00 0.00 2014-08-14 03:54 0.45 12 2014-08-14 03:54 0.37 2014-08-14 16:14 0.51 13 2014-08-13 00:00 0.00 2014-08-13 15:00 0.63 2 4 2014-08-13 15:00 0.13 2014-08-14 07:06 0.67 10 2014-08-14 07:06 1.09 2014-08-14 22:57 0.66 2 3 3 2014-08-13 09:30 0.00 2014-08-13 10:15 0.03 4 7 2014-08-13 07:00 0.00 2014-08-14 07:48 1.03 5 2014-08-13 05:00 0.00 2014-08-13 12:54 0.33 5 6 2014-08-13 20:00 0.00 2014-08-13 23:36 0.15 9 2014-08-13 23:36 0.03 2014-08-14 10:27 0.45 8 2014-08-14 10:27 0.69 2014-08-14 21:41 0.47 6 14 2014-08-13 00:00 0.00 2014-08-13 16:36 0.69 2 2014-08-13 19:00 0.00 2014-08-14 14:15 0.80 1 2014-08-13 01:00 0.00 2014-08-13 15:47 0.62 1 4 2014-08-13 15:47 0.16 2014-08-14 04:40 0.54 10 2014-08-14 04:40 0.99 2014-08-14 17:21 0.53 13 2014-08-13 00:00 0.00 2014-08-13 15:00 0.63 2 6 2014-08-13 20:00 0.00 2014-08-14 00:30 0.19 12 2014-08-14 00:30 0.23 2014-08-14 15:54 0.64 3 3 3 2014-08-13 09:30 0.00 2014-08-13 10:15 0.03 4 4 5 2014-08-13 05:00 0.00 2014-08-14 00:45 0.82 9 2014-08-14 00:45 0.07 2014-08-15 03:51 1.13 7 2014-08-13 07:00 0.00 2014-08-13 16:56 0.41 5 11 2014-08-13 17:00 0.00 2014-08-14 03:54 0.45 8 2014-08-14 03:54 0.41 2014-08-14 15:08 0.47 6 14 2014-08-13 00:00 0.00 2014-08-13 16:36 0.69 2 2014-08-13 19:00 0.00 2014-08-14 14:15 0.80 1 2014-08-13 01:00 0.00 2014-08-13 15:47 0.62 1 12 2014-08-13 19:00 0.00 2014-08-14 07:20 0.51 4 2014-08-14 07:20 0.81 2014-08-14 20:13 0.54 13 2014-08-13 00:00 0.00 2014-08-13 15:00 0.63 2 10 2014-08-13 15:00 0.42 2014-08-14 06:51 0.66 9 2014-08-14 06:51 0.33 2014-08-14 20:24 0.56 4 4 5 2014-08-13 05:00 0.00 2014-08-14 00:45 0.82 7 2014-08-13 07:00 0.00 2014-08-13 16:56 0.41 5 5 8 2014-08-13 18:00 0.00 2014-08-14 05:14 0.47 11 2014-08-14 05:14 0.51 2014-08-14 16:08 0.45 6 2014-08-14 16:08 0.84 2014-08-14 19:44 0.15 14 2014-08-13 00:00 0.00 2014-08-13 16:36 0.69 6 3 2014-08-13 16:36 0.30 2014-08-13 17:15 0.03 2 2014-08-13 19:00 0.00 2014-08-14 14:15 0.80 Note: the unit of the waiting time and operating time is days. 308 Advances in Production Engineering & Management 11(4) 2016 Applying multi-phase particle swarm optimization to solve bulk cargo port scheduling problem 4.3 Discussion On the condition of small samples (Vn = 14, Bn = 6), for Model 1, the total service time, the maximum waiting time and the maximum operating time are minimal, and the makespan is the earliest For Model 2 and Model 3, the total service time, the maximum waiting time and the maximum operating time will be increased slightly, and the makespan can be delayed. For model 4, the total service time is the highest, but the maximum waiting time and the maximum operating time are relatively small. The makespan is also the earliest (Table 11). On the condition of big samples (Vn = 24, Bn = 6), for Model 1, the total service time is still the smallest, but the maximum waiting time and the maximum operating time are increased. Additionally, the makespan is the latest. For Model 2 and Model 3, the total service time is still large, but the maximum waiting time and the maximum operating time are relatively small, and the makespan is competitive compared with Model 1. For model 4, the total service time is still the highest, but the maximum operating time, especially the maximum waiting time, is the lowest. The makespan is the earliest (Table 11). On one hand, the total service time of berth needs to be decreased to satisfy the general requirements, so the model considering priority offers a great advantage. On the other hand, the total service time of berths reflects the capacity; the longer the service time is, the higher the berth service capacity will be. The departure time of the last vessel needs to considered, so the model minimizing makespan offers a great advantage. In short, we can choose a suitable model according to the actual situation in bulk cargo port scheduling. Table 11 Comparative results of four models Model Total service time Maximum waiting time Maximum operating time Makespan Small Big Small Big Small Big Small Big samples samples samples samples samples samples samples samples 1 9.55 24.01 0.56 2.20 0.82 1.04 2014-08-14 20:24 2014-08-16 00:25 2 9.79 24.39 1.09 2.09 1.03 1.09 2014-08-14 22:57 2014-08-15 19:18 3 9.81 24.30 0.99 2.08 0.82 1.09 2014-08-15 03:51 2014-08-15 19:18 4 10.54 27.54 0.84 1.52 0.82 1.09 2014-08-14 20:24 2014-08-15 17:49 5. Conclusion Social and economic progress drives the development of bulk cargo port transportation, so the transportation demand is also growing. However, scheduling management and efficiency are incompatible with its development. The problems are as follows: (1) The scheduling is dispersed. Due to multi-level allocation, vessel scheduling is cumbersome, and much repetitive work is done, thereby affecting the efficiency of scheduling. (2) The management mode lags behind. People basically rely on experience to manage berth scheduling. Therefore, reducing the total service time of berths and maximizing their operational capabilities (minimizing the makespan) are of great importance to improving the efficiency and management of bulk cargo port scheduling. This paper sorts the factors affecting bulk cargo port scheduling, such as the number of vessels, the number of berths, vessel-berthing constraints, the service priority, and the makespan, and establishes the NP model, which aims to minimize the total service time and makespan, and solves the model based on the MPPSO algorithm. Through computational experiments, we obtain information not only about the berth allocation (the starting time, finishing time, waiting time, and operating time) but also the comparative results of different modes (the total service time, maximum waiting time, maximum operating time, and makespan). We can choose a suitable model according to the actual situation in bulk cargo port scheduling. Of course, like all studies, this paper has certain limitations and deficiencies. To guarantee the computability of the model, this paper assumes that the vessels cannot move to other berths if the berth allocation is complete, and the time of vessels moving is ignored. We should further study these topics. Advances in Production Engineering & Management 11(4) 2016 309 Tang, Gong, Liu, Zhang Acknowledgement The study is supported by a project funded by the China Postdoctoral Science Foundation (2016M591194) and the National Natural Science Foundation (71132008,71390334). We appreciate their support very much. References [1] Cheong, C.Y., Tan, K.C., Liu, D.K., Lin, C.J. (2010). Multi-objective and prioritized berth allocation in container ports, Annals of Operations Research November, Vol. 180, No. 1, 63-103, doi: 10.1007/s10479-008-0493-0. [2] Conway, R.W., Maxwell, W.L., Miller, L.W. (2003). Theory of scheduling, Dover publication, New York, USA. [3] Gantt, H. (1919). Organization for Work, Allen and Unwin, London, UK. [4] Xu, D., Li, C.L., Leung, J.Y.-T. (2012). Berth allocation with time-dependent physical limitations on vessels, European Journal of Operational Research, Vol. 216, No. 1, 47-56, doi: 10.1016/j.ejor.2011.07.012. [5] Raa, B., Dullaert, W., Van Schaeren, R. (2011). An enriched model for the integrated berth allocation and quay crane assignment problem, Expert Systems with Applications, Vol. 38, No. 11, 14136-14147, doi: 10.1016/j.eswa. 2011.04.224. [6] Imai, A., Nishimura, E., Papadimitriou, S. (2001). The dynamic berth allocation problem for a container port, Transportation Research PartB: Methodological, Vol. 35, No. 4, 401-417, doi: 10.1016/S0191-2615(99)00057-0. [7] Guan, Y., Xiao, W.Q., Cheung, R.K., Li, C.L. (2002). A multiprocessor task scheduling model for berth allocation: Heuristic and worst-case analysis, Operations Research Letters, Vol. 30, No. 5, 343-350, doi: 10.1016/S0167-6377(02) 00147-5. [8] Arango, C., Cortés, P., Muñuzuri, J., Onieva, L. (2011). Berth allocation planning in Seville inland port by simulation and optimization, Advanced Engineering Informatics, Vol. 25, No. 3, 452-461, doi: 10.1016/j.aei.2011.05.001. [9] Galzina, V., Lujic, R., Saric, T. (2012). Adaptive fuzzy particle swarm optimization for flow-shop scheduling problem, Technical Gazette - Tehnicki vjesnik, Vol. 19, No. 1, 151-157. [10] Meisel, F., Bierwirth, C. (2009). Heuristics for the integration of crane productivity in the berth allocation problem, Transportation Research Part E: Logistics and Transportation Review, Vol. 45, No. 1, 196-209, doi: 10.1016/j.tre.2008.03.001. [11] Golias, M.M., Saharidis, G.K., Boile, M., Theofanis S., Ierapetritou M.G. (2009). The berth allocation problem: Optimizing vessel arrival time, Maritime Economics & Logistics, Vol. 11, No. 4, 358-377, doi: 10.1057/mel.2009.12. [12] Türkogullan, Y.B., Ta§kin, Z.C., Aras, N., Altinel, I.K. (2014). Optimal berth allocation and time-invariant quay crane assignment in container terminals, European Journal of Operational Research, Vol. 235, No. 1, 88-101, doi: 10.1016/j.ejor.2013.10.015. [13] Imai, A., Sun, X., Nishimura, E., Papadimitriou, S. (2005). Berth allocation in a container port: Using a continuous location space approach, Transportation Research Part B: Methodological, Vol. 39, No. 3, 199-221, doi: 10.1016/ j.trb.2004.04.004. [14] Legato, P., Mazza, R.M. (2001). Berth planning and resources optimization at a container terminal via discrete event simulation, European Journal of Operational Research, Vol. 133, No. 3, 537-547, doi: 10.1016/S0377-2217(00)00200-9. [15] Robenek, T., Umang, N., Bierlaire, M., Ropke, S. (2014). A branch-and-price algorithm to solve the integrated berth allocation and yard assignment problem in bulk ports, European Journal of Operational Research, Vol. 235, No. 2, 399-411, doi: 10.1016/j.ejor.2013.08.015. [16] Nishimura, E., Imai, A., Papadimitriou, S. (2001). Berth allocation planning in the public berth system by genetic algorithms, European Journal of Operational Research, Vol. 131, No. 2, 282-292, doi: 10.1016/S0377-2217(00)00128-4. [17] Imai, A., Zhang, J.T., Nishimura, E., Papadimitriou, S. (2007). The berth allocation problem with service time and delay time objectives, Maritime Economics & Logistics, Vol. 9, No. 4, 269-290, doi: 10.1057/palgrave.mel. 9100186. [18] Fu, Y.M., Diabat, A., Tsai, I-T. (2014). A multi-vessel quay crane assignment and scheduling problem: Formulation and heuristic solution approach, Expert Systems with Applications, Vol. 41, No. 15, 6959-6965, doi: 10.1016/j.eswa. 2014.05.002. [19] Wang, J.F., Kang, W.L., Zhao, J.L., Chu, K.Y. (2016). A simulation approach to the process planning problem using a modified particle swarm optimization, Advances in Production Engineering & Management, Vol. 11, No. 2, 77-92, doi: 10.14743/apem2016.2.211. [20] Mocnik, D., Paulic, M., Klancnik, S., Balic, J. (2014) Prediction of dimensional deviation of workpiece using regression, ANN and PSO models in turning operation, Technical Gazette - Tehnicki vjesnik, Vol. 21, No. 1, 55-62. [21] Ting, C.J., Wu, K.C., Chou, H. (2014). Particle swarm optimization algorithm for the berth allocation problem, Expert Systems with Applications, Vol. 41, No. 4 (Part 1), 1543-1550, doi: 10.1016/j.eswa.2013.08.051. 310 Advances in Production Engineering & Management 11(4) 2016 Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 | pp 311-323 http://dx.doi.Org/10.14743/apem2016.4.229 ISSN 1854-6250 Journal home: apem-journal.org Original scientific paper A case-based reasoning approach for non-traditional machining processes selection Boral, S.a, Chakraborty, S.a* department of Production Engineering, Jadavpur University, Kolkata, West Bengal, India A B S T R A C T A R T I C L E I N F O To sustain in the modern era of rapid manufacturing development, it becomes necessary to generate complex shapes on materials which are highly temperature and corrosion resistant, hard to machine, and have high strength-to-weight ratio. Generation of complex shapes on those materials using conventional machining processes ultimately affects surface finish, material removal rate, accuracy, cost, safety etc. Non-traditional machining (NTM) processes have the capability to machine those advanced engineering materials with satisfactory results. But, selection of the most appropriate NTM process for a particular machining application is often a complicated task. Case-based reasoning (CBR), a domain of artificial intelligence, is a paradigm for reasoning new problems from the past experience. In CBR, a memory model is assumed for representing, indexing and organizing past similar cases, and a process model is supposed for retrieving and modifying the past cases and assimilating the new ones. This paper primarily focuses on the application of CBR approach for NTM process selection. Based on different process characteristics and process parameter values, the past similar cases are retrieved and reused to solve a current NTM process selection problem. For this, a software prototype is developed and three real time examples are cited to illustrate the application potentiality of CBR system. © 2016 PEI, University of Maribor. All rights reserved. Keywords: Non-traditional machining processes Process selection Artificial intelligence Case-based reasoning *Corresponding author: s_chakraborty00@yahoo.co.in (Chakraborty, S.) Article history: Received 7 June 2016 Revised 25 October 2016 Accepted 12 November 2016 1. Introduction With the development of newer materials having improved thermal, mechanical and chemical properties, it has now become quite difficult to machine those materials using conventional machining processes. These processes, generally based on cutting and abrasion mechanism, incur higher machining cost while generating complex shape features on composites, ceramics and other advanced engineering materials. The achieved surface quality and dimensional accuracy of the machined components are also not satisfactory, and often fail to meet the desired target In these machining processes, unwanted material from the parent workpiece is generally removed employing mechanical energy. This energy is supplied by means of a cutting tool kept in contact with the workpiece, causing shear deformation along the shear plane, leading to chip formation. New exotic work materials and complex geometrical shapes on those materials have been putting more pressure on the capabilities of the conventional machining processes. This leads to the development and deployment of a new set of machining processes, popularly known as non-traditional machining (NTM) processes. In these processes, unwanted material is removed from the parent workpiece using various forms of energy, like chemical, thermal, mechanical, electrical or combination of those energies. In an NTM process, there is no direct contact between the 311 Boral, Chakraborty cutting tool and the workpiece. In abrasive jet machining process, excess material is removed by means of microscopic chips and in electrochemical machining process by electrolytic dissolution. In laser beam machining process, there is even no need of any cutting tool. It is also not necessary that the cutting tool should be harder than the workpiece material in an NTM process. Now-a-days, it has become easier to generate complex shapes on materials, like steel, carbide, titanium and its alloys, ceramics, superalloys (Inconel 718, hastelloy) etc. employing NTM processes [1,2]. Till date, there have been approximately 20 NTM processes developed and applied in modern manufacturing industries. Selection of the best suitable NTM process for a particular work material and shape feature combination is generally made by a domain expert on the basis of various factors, such as workpiece material, shape feature to be generated, material removal rate, surface finish, surface damage, corner radii, tolerance, cost, safety, power requirement etc. Thus, an expert in this domain must have a vast and in-depth knowledge about the characteristics and capabilities of different available NTM processes. But, in the present manufacturing scenario, most of the process engineers lack the requisite domain knowledge and availability of experts is also sometimes constrained. Usually, a domain expert acquires knowledge from the past experience as well as from other reliable sources. Taking this concept as a plinth, when an expert attempts to select an NTM process for a given machining application, he/she just recalls the similar past situations and their solutions. Thus, based on the similar past problems and their solutions, new NTM process selection cases are solved. This entire cognitive process of a domain expert's thinking has given birth to a new branch of artificial intelligence (AI) technique, known as case-based reasoning (CBR) approach. This CBR approach is applied here for NTM process selection. In this paper, in order to choose the most suitable NTM process for a specific machining application, an exhaustive case-base containing the machining characteristics of various available NTM processes and their pertinent process parameters is first created. These machining characteristics and process parameter data are later used to select the feasible NTM processes according to the end requirements. The selection procedure is based on retrieval of the best matched case from the case-base using the nearest neighbourhood technique, while calculating the similarity score between two cases. The best matched case, which is retrieved from the case-base according to the values of different process characteristics as set by the process engineer/end user, has the similarity score greater than the other cases. To automate and simplify the application of CBR approach in NTM process selection, a software prototype having a graphical user interface (GUI) is designed and developed in Visual Basic 6.0. The developed system simultaneously considers both the user requirements (product characteristics) and technical requirements (process characteristics) for a given NTM process selection problem. 2. Literature review Using two multi-attribute decision making (MADM) tools, i.e. analytic hierarchy process (AHP) and technique for order preference by similarity to ideal solution (TOPSIS), Yurdakul and Cfogun [3] attempted to simplify the NTM process selection procedure for the manufacturing personnel. A list of feasible NTM processes satisfying the users' requirements was first generated and those processes were then ranked based on their suitability to meet the desired machining operation. An expert system was developed by Chakraborty and Dey [4] for selecting the best NTM process under constrained material and machining conditions. It would rely on the priority values of different criteria and sub-criteria for a specific NTM process selection problem, and the NTM process with the highest acceptability index was finally identified. Chakraborty and Dey [5] developed a quality function deployment (QFD)-based expert system for NTM process selection. An overall score for each of the NTM processes was estimated using the weights extracted from the house of quality matrix for various process characteristics. The overall scores of some of the NTM processes simultaneously satisfying certain critical criteria requirements were again compared and the NTM process having the maximum score was finally selected as the optimal choice. A web-based knowledge base system was proposed by Edison Chandraseelan et al. [6] for identifying the most suitable NTM process to meet some input parametric requirements, 312 Advances in Production Engineering & Management 11(4) 2016 A case-based reasoning approach for non-traditional machining processes selection like type of the work material, shape application, process economy, and other process capabilities, like surface finish, corner radii, width of cut, tolerance etc. Sadhu and Chakraborty [7] applied an input minimized Charnes, Cooper and Rhodes (CCR) model of data envelopment analysis for NTM process selection. Employing weighted-overall efficiency ranking method of MADM theory, the efficient NTM processes were ultimately ranked in descending order of their priorities. Temufin et al. [8] designed a fuzzy decision support model for NTM process selection while assessing the potentials of some distinct NTM processes. Chatterjee and Chakraborty [9] proved the application potentiality of evaluation of mixed data (EVAMIX) method for solving NTM process selection problems using three demonstrative examples. Roy et al. [10] integrated fuzzy AHP and QFD techniques for selection of NTM processes based on some predefined customers' perspectives. Temufin et al. [11] solved the NTM process selection problem under fuzzy and crisp environment, and proposed a decision support model to guide the process engineers to explore the potentials of some distinct NTM processes. The applicability of the proposed model was also validated. Khandekar and Chakraborty [12] applied fuzzy axiomatic design principles for selection of NTM processes. Madic et al. [13] demonstrated the applicability, suitability and computational procedure of operational competitiveness ratings analysis (OCRA) method for solving NTM process selection problems. Nowadays, CBR as a part of cognitive science, has been emerged out as an interesting research topic. Amen and Vomacka [14] employed CBR approach as a tool for selection of material and heat treatment process from an exhaustive database to simplify the task of a designer. Khe-mani et al. [15] applied CBR approach in fused cast refractory manufacturing industry. Fang and Wong [16] applied a hybrid CBR approach in agent-based negotiation for effective supply chain management Armaghan and Renaud [17] adopted CBR approach to prove the complementary nature of multi-criteria decisions and CBR approach. Although the past researchers applied numerous MADM methods and developed different distinct decision aids for selection of NTM processes for varying machining applications, but till date, no attempt has been put forward on selection of NTM processes using CBR approach. This paper thus proposes development of a decision making model based on CBR approach for selecting the best suited NTM process for a given machining application. It is observed that CBR is the correct and simplest approach in this domain where availability of experts is sometimes constrained. In CBR approach, a set of feasible NTM processes is first retrieved from the case-base satisfying the work material and shape feature combination. Based on the user and technical requirements, it then identifies the best matched NTM process from the stored similar cases. The past cases are just reused here for NTM process selection for providing the optimal solution. 3. CBR approach Intelligence, being a part of cognitive science, can be defined as the process involving rational and abstract thinking. It is often goal oriented and purposeful. It consists of knowledge and feats, both conscious and unconscious, which are acquired through continuous study and experience. The AI is actually the intelligence in machines. Intelligent system is the basement of knowledge engineering. It involves several tasks, like knowledge acquisition, creation of a knowledge base, knowledge representation and use of the acquired knowledge. The represented knowledge is basically used for reasoning or inference. In AI, knowledge is represented using symbols along with heuristics or rules of thumb. While using these heuristics, one should not have to rethink when a similar problem is encountered. The expert system can be defined as an intelligent computer program that uses knowledge and inference procedures to solve problems that are difficult enough requiring significant human expertise for their solution. Basically in expert system, knowledge is represented using 'if-then' rules. The CBR approach is a part of AI technique that utilizes information stored in the knowledge base, when similar past problems are encountered again. It provides solution to the present problem that is almost similar to the past. In CBR approach, a problem is represented as an input in the present situation. It just retrieves the most similar case to the new one from its case-base while calculating the similarity score over the defined parameters. It first searches the case history Advances in Production Engineering & Management 11(4) 2016 313 Boral, Chakraborty and chooses that case having the closest similarity to the current problem. In CBR system, the case-base is well structured and documented. The case representation may be flat, where all cases are represented at the same level, or it can be hierarchical, expressing relationship between cases and sub-cases. There are four major steps that constitute a CBR system, i.e. retrieve, reuse, revise and retain. Thus, it is also called as 4-R cycle or CBR cycle, as shown in Fig. 1. When a problem occurs in the current situation, similar past situations are retrieved from the case-base. Reusing the past cases, a predictable solution to the current problem is thus provided. If there is a need of any revision, the retrieved data are revised and retained as a new case in the case-base for future use [18-21]. New case -------------- Confirmed [Retain / ) Retrieve Most solution " * Case-Base -► similar cases Fig. 1 A CBR cycle or 4-R cycle Retrieving the most similar case along with the solution is based on some logical expressions. The similarity between two cases is usually measured with respect to each parameter. It also depends on the type of parameter (beneficial or non-beneficial) being used. The followings are the most common methods for calculating similarity between two cases: a) Numeric: Sim (a,b) = |a - bl/Range where Range is the difference between the upper and lower boundaries of a set. b) Symbolic: Sim(a,b) = 1 if a = b = 0 if a * b (1) c) Multi-valued: Sim(a, b) = Card (a) n Card (b) Card (a) U Card (b) where Card is the cardinality (size) of a set. d) Taxonomy: (2) (3) Sim(a, b) = h (common node (a, b)) (4) min(h(a),h(b)) where h is the height (number of levels) of the specified taxonomy tree. The procedural steps of a CBR approach are presented as below: a) A solution is first defined using several parameters. One of the parameters should be chosen carefully so that it would remain unique throughout the documentation procedure, e.g. case number. b) A huge set of known solutions is put into the case-base of CBR system. An existing database can also be used for this purpose. c) The CBR system generally reads the database and organizes a copy of its own. 314 Advances in Production Engineering & Management 11(4) 2016 A case-based reasoning approach for non-traditional machining processes selection d) The user generally formulates a query according to the end requirements. All the available variables are first displayed. The user has the option to choose all or few variables based on the problem statement. The query includes those variables as set by the user. The user also has the option to allocate different priority weights to the considered variables. e) As a result of the user-defined query, CBR system may display a number of cases or the best matched case. It may also be possible that none of the cases would match the query exactly. Favouring CBR technique as the most efficient tool for NTM process selection is a challenging task, as several other approaches have already been available for the same purpose. It is observed from the available literature that none of the MADM methods, like AHP, EVAMIX, TOPSIS etc. can provide complete solution when the domain is ill-structured and murky. The working principle of CBR is based on some available specific experiences instead of abstracted rules. It is considered as a useful tool if the utilization of prior experience is more vital than to produce a thoroughly optimized solution according to the specifications. The CBR approach has no optimizing potentiality, but it can be used for searching, not for calculations. Its efficiency is determined by fast retrieval of the most similar cases from the case-base. The principle of CBR also states that it can find the similarities between cases but not reasons. So, it is unable to judge how important the encountered departures are that can be determined only by an experienced user. A comparison between the existing search techniques and the adopted CBR approach is elucidated in Table 1. Table 1 Comparison between different search techniques and CBR approach Method Flexibility Operational approach Computational time Programming complexity Decision maker's involvement Type of data Genetic algorithm Medium (lack of learning ability) Artificial neural High network Simulated annealing Expert system Medium Medium CBR High Population based probabilistic search and optimization technique High (based on the desired accuracy and termination criterion) Mimics the working principle of biological High neurons Cooling process of molten metal is modeled artificially to construct an optimization algorithm Exact matching of input and stored data producing several Medium 'if-then' rules for inference Notion of similarity between present and prior stored cases Low High Medium (based on the desired accuracy and termination criterion) Medium High Medium Low High High High Medium Numerical Medium Numerical Numerical Both numerical and textual Both numerical and textual 4. CBR-based approach for NTM processes selection Although CBR approach has already been successfully applied in various fields of mechanical engineering, such as material selection, design selection, parts selection for automobile industries etc., no attempt has still been made for its application in the domain of NTM process selection. The CBR approach has the potential to provide complete information about a case where minimum information is available to the user. It yields the best results when the user provides detailed query information. Advances in Production Engineering & Management 11(4) 2016 315 Boral, Chakraborty While selecting the most suitable NTM process for a particular machining application, the process engineer has to consider several machining characteristics of the available NTM processes. In the developed CBR approach-based decision making model, nine NTM processes, i.e. abrasive jet machining (AJM), abrasive water jet machining (AWJM), electric discharge machining (EDM), laser beam machining (LBM), ultrasonic machining (USM), electrochemical machining (ECM), electrochemical discharge machining (ECDM), plasma arc machining (PAM) and wire electric discharge machining (WEDM) are taken into consideration. As the process characteristics, type of the workpiece material, shape feature to be generated, material removal rate (MRR) (in mg/min), surface roughness (SR) (in [im), surface damage (SD) (in [im), tolerance (Tol) (in mm), overcut (OC) (in mm), corner radii (CR) (in mm), taper (TP) (in mm/mm), cost (C) (in relative (R) priority scale), power (P) (in kW) and safety (S) (in R scale) are considered. For cost, the R scale is set as 1 - lowest, 2 - very low, 3 - low, 4 - medium, 5 - high, 6 - very high and 7 - highest. On the other hand, for safety, the R scale is set as 1 - highly safe, 2 - safe and 3 - attention required. As work materials, a) aluminium, b) aluminium alloys, c) ceramics, d) composites, e) glass, f) steel, g) superalloys and h) titanium are considered in this model. The above-mentioned NTM processes can generate a) hole (precision) (0.03 mm < D < 0.13 mm), b) hole (standard) (L/D < 20), c) hole (standard) (L/D > 20), d) through cut (shallow) (t/w < 2), e) through cut (deep) (t/w > 2), f) through cavity (standard) (t/w > 10), g) through cavity (precision) (t/w < 10), h) pocket (shallow) (t < 1 mm), i) pocket (deep) (t > 1mm) and j) surface of revolution feature on the work material (where L is the length of the hole, D is the diameter of the hole, t is the thickness and w is the width of the machined feature). The relevant machining characteristics data for different NTM processes are accumulated from experimentations, machining data handbooks and other reliable resources to create the corresponding case-base. The collected data are then organized in a structured manner in MS Access. The step-wise operational procedures of the developed CBR system for selecting the best suited NTM process for a particular machining application are depicted as follows: Step 1: When the developed CBR system starts to run, the first screen, as shown in Fig. 2, appears to the end user where the type of the work material to be machined and type of the shape feature to be generated can be chosen from the given options as the primary inputs to the system. Step 2: After clicking 'OK' button, a list of feasible NTM process(es) capable of generating the desired shape on the specified work material is displayed. For this, Eqn. (2) is utilized for filtering and retrieving the data. Step 3: When the user presses 'Next' button, another screen, as shown in Fig. 3, is displayed to identify the most suitable NTM process from the list of feasible processes while satisfying the set machining requirements. Step 4: In this screen, the end user has to choose the desired process characteristics based on which the final NTM process selection is made. Step 5: When 'Enter range' functional key is clicked, the required number of empty cells are automatically generated where the input ranges for the selected NTM process characteristics can be provided. Step 6: After inputting the desired ranges of values, pressing of 'Best NTM process' button identifies the most suitable NTM process for the specified machining application while satisfying the set criteria values. For retrieving the best NTM process in this step, Eqn. (1) is employed. Step 7: The actual retrieved values of all the technical characteristics for the best matched NTM process are also displayed. Step 8: When 'Best NTM process' button is clicked, the technical details (tentative settings of the associated process parameters) of the best matched NTM process are also available, as shown in Fig. 4. 316 Advances in Production Engineering & Management 11(4) 2016 A case-based reasoning approach for non-traditional machining processes selection Although in step 5, there is an option for entering ranges of process characteristic values, but if the developed CBR system does not find any data within those ranges from the case-base, it would retrieve the possible data nearest to the query set. For a particular NTM process selection problem, MRR is the sole beneficial attribute where its value is always required to be maximized. On the other hand, SR, SD, Tol, TP, OC, CR, C, P and S are non-beneficial attributes requiring their minimum values. The best matched case should have the highest similarity score, which is calculated with respect to each of the process characteristics. After summing up these similarity scores for the set process characteristics for each case, the NTM process having the highest similarity score is chosen as the most suitable option. 5. Illustrative examples 5.1 Example 1: Standard hole on composite material In this example, standard holes are to be generated on a composite material. After providing the inputs of composite as the work material and hole (standard) as the shape feature options in the primary selection window of Fig. 2, a set of feasible NTM processes consisting of AJM, AWJM, ECDM, ECM, EDM, LBM and USM is displayed, when 'OK' button is clicked. All the processes can generate standard holes on composite materials. In the next window of Fig. 3, MRR, SR, Tol, OC, CR and C are opted as the most important process characteristics based on which the final NTM process selection is to be made. In this example, the desired input ranges for those process characteristics are set as MRR 100-1000 mg/min, SR 2-12 ^m, Tol 0-0.5 mm, OC 0-0.05 mm, CR 0-0.5 mm and C 1-4 (in R scale). Now, when 'Best NTM process' functional button is clicked, LBM process is identified as the best matched case, capable of meeting the set process characteristic values. It is interesting to observe that apart from the set process characteristics, values of the other process characteristics are also available for the best matched NTM process. In this example, the selected LBM process can achieve values of MRR as 286.08 mg/min, SR as 2.63 [j.m, SD as 102 [j.m, Tol as 0.02 mm, OC as 0.001 mm, CR as 0.5 mm, TP as 0.05 mm/mm, C as 1 (in R scale), P as 0.23 kW and S as 3 (in R scale). In Fig. 4, the process engineer can also have an idea about the settings of different machining parameters of LBM process. These are the tentative process parametric settings and for achieving the maximum machining performance, fine-tuning of these settings is often necessary. A real time photograph of LBM process is also available in Fig. 4. Fig. 2 Primary selection window for Example 1 Advances in Production Engineering & Management 11(4) 2016 317 Boral, Chakraborty Final selection window using CBR approach Process characteristics w Material removal rate (MRR) P Surface roughness (SR) r Surface damage (SD) I? Tolerance (Tol) P Overcut (OC) P Corner radii (CR) r Taper (TP) P Cost (C) r Power (P) r Safety (S) Select process characteristics and input ranges Material removal rate (mg/min) Surface roughness (pm) Tolerance (mm) Overcut (mm) Corner radii (mm) Cost (1-7) (R scale) Enter range 0.05 Best NTM process LBM R Scale for cost : 1 ^Lowest, 2- Very low, 3= low, 4= Medium, 5 = High, 6= Very high, 7= Highest R Scale for safety : 1 -Highly safe. 2- Safe. 3- Attention required Matched case : Process MRR SR SD Tol OC CR TP C P S 1 LBM 1286.08 2.63 1 102 0.02 0.001 1 0.Q5 1 0.05 I I 1 I 0.23 I 3 I Fig. 3 Best NTM process for Example 1 0 LBM El Laser beam machining (LIBM1 Pulse power (kW) Gas pressure (Bar) Cutting speed (mm/min) Stand off distance (mm) 0.23 9.8 0.5 Pulse frequency (Hz) Pulse width (ms) Nozzle dia (mm) Focal length (mm) i ^ " Process |_BM Workpiece material Composite isl si Exact material AA7075/SÍC MMC m * ^ Gas used Nitrogen Shape Hole (Standard) [LiD<20] 230 0.2 1.2 116 Fig. 4 Details of LBM process 5.2 Example 2: Standard through cavity on ceramics Here, the process engineer wants to generate a standard through cavity on a ceramic work material. In the primary selection window, as shown in Fig. 5, the developed CBR approach first extracts five NTM processes, i.e. AJM, AWJM, EDM, USM and WEDM as the feasible options satisfying the said work material and shape feature combination requirement. In Fig. 6, MRR, SR, Tol, OC, CR, C and S are chosen as the most important process characteristics based on which the final NTM process needs to be selected. Based on the ranges of values for these process characteristics, USM process is identified as the best matched case for this machining application. For USM process, the attainable process characteristics are MRR as 131.96 mg/min, SR as 0.66 [im, SD as 25 [im, Tol as 0.014 mm, OC as 0.15 mm, CR as 0.08 mm, TP as 0.005 mm/mm, C as 5 (in R scale), P as 0.4 kW and S as 1 (in R scale). 318 Advances in Production Engineering & Management 11(4) 2016 A case-based reasoning approach for non-traditional machining processes selection El. Primary selection window ^ i ■Work material- Ceramic Shape Through cavity (Standard) [t/w>10] OK Feasible NTM process(es) AJM AW J M EDM USM WEDM Next Fig. 5 Primary selection window for Example 2 In Fig. 7, the tentative parametric settings and the technical specifications of USM process along with its actual photograph are displayed to guide the process engineer to achieve the best machining performance. 5 Final selection wl ndow us i ng CBRapp roach H I i Select process characteristics and input ranges Process characteristics P Material removal rate (MRR) w Surface roughness (SR) r Surface damage (SD) P Tolerance (Tol) W Overcut (OC) F Corner radii (CR) r Taper (TP) W Cost (C) r Power (P) W Safety (S) Enter range Material removal rate (mg/min) 10 100 Surface roughness (pm) 2 20 Tolerance (mm) 0 0.05 Overcut (mm) 0 0.01 Corner radii (mm) 0 0.05 Cost (1-7) (R scale) 1 4 Safety (1-3) (R scale) 1 2 i Best NTM process! USM R Scale for cost : 1=Lowest, 2- Very low, 3= low, A- Medium, 5 = High, 6= Very high, 1- Highest R Scale for safety : 1=Highly safe, 2- Safe, 3= Attention required Matched case : Process MRR SR SD Toi OC CR TP USM 131.946 0.66 25 0.014 0.15 0.08 0.005 0.4 Fig. 6 Best NTM process for Example 2 Advances in Production Engineering & Management 11(4) 2016 319 Boral, Chakraborty a usM Ultrasonic machining (USM) Process Shape |USM Through cavity (Standard) [t/w>10] Workpiece material Exact material Tool material Abrasive used Slurry media Feed rate (mm/min) Amplitude of vibration (|jm) Ceramic Zirconia Stainless Steel Boron Carbide Water 1.08 32 Slurry concentration (%) Frequency of vibration (Hz) Material thickness (mm) 20 Fig. 7 Details of USM process 5.3 Example 3: Shallow through cutting on steel In this example, a shallow through cutting operation needs to be performed on a standard steel plate. For this work material and shape feature combination, the CBR system first recognizes AJM, AWJM, ECM, EDM, LBM and PAM as the six feasible NTM processes, as shown in Fig. 8. Then, in Fig. 9, seven process characteristics, i.e. MRR, SR, SD, Tol, OC, CR and C are identified by the process engineer for the final selection of the most suited NTM process for the considered machining application. In this window, the ranges of values of the set process characteristics are also provided. 5 Primary selection window (Si \w » I rWork material- Steel Shape Through cutting (Shallow) [t/w<2] OK Feasible NTM process(es) AJM - AWJM ECM — EDM - LBM - Next Fig. 8 Primary selection window for Example 3 320 Advances in Production Engineering & Management 11(4) 2016 A case-based reasoning approach for non-traditional machining processes selection Final selection window using CBR approach Process characteristics W Material removal rate (MRR) W Surface roughness (SR) R Surface damage (SD) F Tolerance (Tol) W Overcut (OC) W Corner radii (CR) r Taper (TP) I? Cost (C) r Power (P) r Safety 0 or (1 - x/ - A/P1,/) > 0. It is further assumed that all defective items can be reworked and repaired. The rework processes starts immediately after the end of regular production processes in each production cycle (see Fig. 2), at a rate of P2,/. Fig. 1 Inventory level of perfect quality common intermediate parts and customized final products in the proposed two-stage multi-product vendor-buyer integrated inventory system with rework Advances in Production Engineering & Management 11(4) 2016 335 Chiu, Kuo, Chiu, Hsieh t V-T1 *1,1 ll\ tx l rj-t- 'Ti'uiie1 Time Stage 2 Fig. 3 Inventory level of common intermediate parts waiting to be fabricated into customized final products in the stage 2 of the proposed two-stage multi-product system Upon completion of the production in stage 1, L different lots of common intermediate parts are made ready for the production in stage 2. They are fabricated in sequence into customized end products under the common production cycle time policy. The inventory level of common intermediate parts waiting to be fabricated in stage 2 is depicted in Figure 3. In stage 2, after the completion of rework process (t2,i) of each end product i, fixed quantity n installments of the finished batch are transported to customers at a fixed interval of time in the delivery time t3,i (see Fig. 1). The inventory level of end products at the buyers' side during a production cycle is depicted in Figure 4 (which is similar to Fig. 3 in [9]). The following are additional notation used in this study (where i = 1, 2,..,L, represents L different products in stage 2; and i = 0 denotes the common intermediate part in stage 1): T - Production cycle length, one of the decision variables, n - Number of fixed quantity installments of the finished batch to be delivered in each cycle, the other decision variable, a - Completion rate of common intermediate part as compared to the finished product, Qi - Production lot size for product i, Ki - Production setup cost for product i in a production cycle, Ci - Unit production cost for product i, hu - Unit holding cost for product i, h2i - Holding cost per reworked item for product i, D{ 11 D D, m Time t i tv1 H ' ui z<('l,l l'zi¡_ ■ I. .1. . I. t '/ t 1,0 i 2.0 i,i n 21 A-/, =x,t «¡V XJtt^ 1 Time f l/"L¿.¿L 2| I "T" Time iW^-î __L----j--------^ Fig. 4 Inventory level of customized final products at the buyers' side during a production cycle [9] 336 Advances in Production Engineering & Management 11(4) 2016 Effect of delayed differentiation on a multi-product vendor-buyer integrated inventory system with rework h3,i - Unit holding cost for stocks stored at customer's side, h4,/ - Unit holding cost for safety stocks stored at producer's side, Crj - Unit reworking cost for product i, ty - Production uptime for product i in a production cycle, t2,i - The reworking time for product i in a production cycle, t3i - Delivery time for product i in a production cycle, Hi - Inventory level of common intermediate part at the time of producing end product i, Hi,i - Maximal level of perfect quality items i in the end of regular production, H2,i - Maximal level of perfect quality items i in the end of rework process before delivery, K^ - Fixed delivery cost per shipment for product i, Ctj - Unit delivery cost for product i, tn,i - A fixed interval of time between each installment of finished items of product i to be delivered to customer during downtime t3i, I(t)i - On-hand inventory level of perfect quality items i at time t, Id(t)i - On-hand inventory level of defective items i at time t, Ic(t)I - On-hand inventory level of finished product i at time t, at customer's side, Ii - The left-over number of finished items of product i in each tni, at customer's side, Di - Number of finished items of product i to be transported to customer in each shipment, TC(T, n) - Total production-inventory-delivery cost per cycle, E[T] - The expected production cycle length, E[TC(T, n)] - The expected production-inventory-delivery cost per cycle, E[TCU(T, n)] - The long-run average costs per unit time for the proposed model. 2.1 Modeling and analysis A two-stage EPQ-based production plan considering the postponement is proposed to satisfy annual demand Ai of L different customized products. From Figure 1, we observe the production cycle time as T = tu + t2ii + t3ii for i = 0,1,2,-, L (1) In stage 1, the production lot-size of common intermediate parts Qo, depends on the sum of production lot sizes Qi of L different products to be made in the stage 2. Therefore, we obtain the following equations (refer to Fig. 1): Qi = W for i = 1,2,-,L (2) Qo = Zi=i& = ¿oT (3) ^ =lT- = -T^T (4) ^1,0 ^1,0 "1,0 Wi,o = ti,o(^i,o-di,o) (5) H2,o = HliQ + P2,ot2,o = T.i=iQi (6) £ _ x0 Qo _ ^1,0^1,0 _ ^2,0~^1,0 (7) , ^2,0 ^2,0 ^2,0 tfi = H2,o -QI (8) Hi = H(i.i)-Qi for j = 2,3,-,L (9) Hl = %-i) ~QL = 0 (10) In stage 2, for fabrication of L different products we obtain the following equations directly from Figs. 1 to 4 (where i = 1, 2,..., L): ti,i = ~jT~.= (11) ^ 1,1 r !,j a1L (12) (13) (14) Í? II {Pu H2,Í = --Hi,í + P2,ít2,í _ xíQÍ _ di,í P2.1 P 2.1 P2.1 Advances in Production Engineering & Management 11(4) 2016 337 Chiu, Kuo, Chiu, Hsieh A = — ' n h Di ^ítn,í ni i = Äi(tu + t2,i) (15) (16) (17) (18) 2.2 Cost analysis Inventory holding costs for common intermediate parts (including perfect and imperfect quality items) during ti,o, t2,o, and t3,o, are (see Figs. 1 and 2) hi,o pf^ + ("2,0+^0)t2,° + H=1 Ht(t14 + t2í)] + hli0 [(dl,ot;,o)tl,° (19) In stage 2, inventory holding cost for common intermediate parts waiting to be fabricated into customized end products (see Fig. 3) is (20) Inventory holding costs for imperfect quality items waiting to be reworked in both stages are h2,0 pf^o)] + [h2ji (tw)] (21) In stage 2, fixed and variable delivery costs and inventory holding cost for finished product i waiting to be distributed in t3,i are %=1[nKu + CT>iQi]+ Ef=1 [hxi H2iit3ii} The stock holding cost for end product i stored at customers' sides (see Fig. 4) is XU K + ^htn, + ^f^]} (22) (23) The overall cost per cycle TC(T, n) for the proposed system, includes production setup cost, variable production cost, reworking cost, holding cost, and safety stock cost in both stages; and fixed and variable delivery costs and holding costs for stocks stored at customers' side in stage 2. Hence, TC(T, n) is TC(T,n) = K0 + C0Q0 + Cr qxqQq + h2,o (^m) (t2,o) + h4,0(xQQQ)T +hxo pf^ + ^^ (í2,0)+ ^ (ii,o)+ H=1 Ht(t14 + t2i)] Ki + CiQi + CR,iXiQi+nK1,i+CT,iQi + h2,i(^^^^)(t2,i) " + (24) Substituting Eqs. 1 to 18 in Eq. 24 and taking randomness of defective rate into account, and the long-run average system costs £"[TCU(T, n)] can be derived as follows: E[TCU(T,n)] = E|~TC(T,n)] E[T] + C0À0 + CR0À0E[X0] +Z0T} + + CiAi + Cr,ÍAÍE[xí] + ^i + CT,¿A¿] + h-^{ô2A - 5-f} , | h2}iT^jE[xl]2 h3}lT*.j 2P2,i 2 P1¡ P2i n\ 4,1 1 u (25) where 338 266 Advances in Production Engineering & Management 11(4) 2016 Effect of delayed differentiation on a multi-product vendor-buyer integrated inventory system with rework Zo = 1 | 2E[x0] E[x0]2j | h2>0A20E[x0]2 2,0 2P- 2,0 U L Í=1 7=1 n Su = T.-TT.- ;and S2, = If+ + + h40À0E[x0] E[Xíi Ai 2,1 1,1 2,1 /or i = 1,2,...,L (26) 3. Convexity and the optimal decision Upon obtaining the long-run average system costs E[TCU(T, n)], we then prove it is a convex function by applying the Hessian matrix equations [45] to verify that Eq. 27 holds. [T n] From Eq. 25 we obtain 'd2E[TCU(T,n)] á2g[rCU(r,n)]\ dT2 dTdn ! • |Tl d2E[TCU(T,n)] d2E[TCU(T,n)] I LnJ >0 dTdn dn2 dT 2P 2,1 1,1 2,1 dE[TCU(T,n)] _ 2K0 dT2 dE[TCU(T,n)] dn dE[TCU(T,n)] dn2 d2E[TCU(T,n)] dTdn Substituting Eqs. 29, 31, and 32 in Eq. 27, we obtain /d2E[TCU(T,n)] d2E[TCU{T,nXP [T n] • = 2Ko + yL (2Ki , 2ngu| = 1 I j,3 J dT2 dTdn d2E[TCU(T,n)] d2E[TCU(T,n)] (27) (28) (29) (30) (31) (32) (33) dTdn dn2 Because K0, K, and T are all positive, we find Eq. 33 is positive. Hence, E[TCU(T, n)] is a strictly convex function for all T and n different from zero. In order to simultaneously determine production-shipment decision for the proposed system, we can solve the linear system of first derivatives of E[TCU(T, n)] with respect to T and n, respectively, by setting these partial derivatives equal to zero. With further derivations we find r = N K0+Zf=1(Kj+nKu) Zo+Sti —{Ô2l-—\+ 2P2j1 JM., +h4,iÄiE[xi] (34) and I2 Zo+Stl- 1 2 K,í P2,i n \ J • > (35) n* = Advances in Production Engineering & Management 11(4) 2016 339 Chiu, Kuo, Chiu, Hsieh 4. Numerical example and discussion The following numerical example is used to show the practical uses of research results obtained in the previous section. Consider a manufacturer must fabricate five different products and they share a common intermediate part that has completion rate a = 0.5 (i.e., halfway done). To ease comparison efforts for readers, we reconsider a numerical example used in a prior study [9] regarding optimization of a single-stage multi-product system without adopting postponement in its production. Annual production rates of five end products Pu = 58,000, 59,000, 60,000, 61,000, and 62,000 units, respectively; annual demands h = 3,000, 3,200, 3,400, 3,600, and 3,800 units, respectively; annual reworking rates PZi = 46,400, 47,200, 48,000, 48,800, and 49,600 units, respectively; setup costs K = $17,000, $17,500, $18,000, $18,500, and $19,000, respectively; unit fabrication costs Ci = $80, $90, $100, $110, and $120, respectively; the defective rates Xi follow uniform distribution over the intervals [0, 0.05], [0, 0.10], [0, 0.15], [0, 0.20], and [0, 0.25], respectively; and unit reworking costs CR,i = $50, $55, $60, $65, and $70, respectively. Based on common intermediate part's completion rate a = 0.5, a straightforward relationship 1/a is assumed for its relevant production rates. Hence, in the proposed two-stage single-machine production scheme we have Pt,0 = (1/a)*(the mean of P1,/s) = 120,000 and P2,0 = (1/a)*(the mean of P2,/s) = 96,000. The relationship between common intermediate part's relevant costs and its completion rate a can either be linear or nonlinear. Both cases are investigated in the following subsections. 4.1 Case 1: Analysis of linear relationship of cost relevant variables If the relationship between practical fabrication related cost of common intermediate part (or called 'the value' of common part) and its completion rate a is linear, then for a = 0.5 we have the following linear-based relevant values of variables in our proposed system: C0 - $40, unit fabrication cost for common intermediate part, K0 - $8,500, setup cost for common intermediate part, Cr,0 - $25, unit reworking cost for common intermediate part, ^1,0 - $5, unit holding cost for common intermediate part, ^4,0 - $5, unit safety stock cost for common intermediate part, h2,0 - $15, unit holding cost for common intermediate part during the reworking processes, K - Setup costs of end products are $8,500, $9,000, $9,500, $10,000, and $10,500 respectively, X0 - [0, 0.04], the interval uniformly distributed defective rate in the production of common intermediate part, Ci - Unit production costs of end products are $40, $50, $60, $70, and $80, respectively, hu - Unit holding costs of end products are $10, $15, $20, $25, and $30, respectively, Pu - Annual production rates of five end products are 112,258, 116,066, 120,000, 124,068, and 128,276 units, respectively; they are simply calculated by Pu = 1/(1/Pu - 1/Pi,0), Xi - End items' defective rates follow the uniform distribution over the intervals [0, 0.01], [0, 0.06], [0, 0.11], [0, 0.16], and [0, 0.21], respectively, CR,i - Unit reworking costs of end products are $25, $30, $35, $40, and $45, respectively, P2,i - Annual reworking rates of five end products are 89,806, 92,852, 96,000, 99,254, and 102,621 units, respectively; they are simply calculated by P2,i = 1/(1/P2,i - 1/P2,0), h2,i - Unit holding cost per reworked items of end products are $30, $35, $40, $45, and $50, respectively, Ku - Fixed delivery costs per shipment: $1,800, $1,900, $2,000, $2,100, and $2,200, respectively, Cr,i - Unit delivery costs of end items are $0.1, $0.2, $0.3, $0.4, and $0.5, respectively, h3j - Unit holding costs at the customer's side are $70, $75, $80, $85, and $90, respectively, h4,0 - Unit safety stock costs of end products are $10, $15, $20, $25, and $30, respectively. 340 Advances in Production Engineering & Management 11(4) 2016 Effect of delayed differentiation on a multi-product vendor-buyer integrated inventory system with rework First, the annual demand for common intermediate parts Xo = 17,000 can be obtained by applying Eqs. 2 and 3. Then, by calculating Eqs. 34, 35, and 25, we derive the optimal number of deliveries n* = 3, optimal production cycle time T* = 0.4614 (years), and the expected system costs per unit time E[TCU(T*, n*)] = $2,145,834. Figure 5 depicts the effects of variations of the production cycle time Ton the expected system costs E[TCU(T, n)]. The behavior of E[TCU(T, n)] with respect to the common intermediate part's completion rate a is exhibited in Figure 6. It can be seen that as the completion rate a increases, the long- run expected system costs E[TCU(T, n)] decreases, and the proposed model realizes a system cost savings of 3.76 % at a = 0.5 (i.e., system costs decreased from $2,229,658 [9] to $2,145,834) as compared to that in prior study which used a single-stage production scheme. This analytical result demonstrates that the proposed two-stage multi-item production scheme with delayed differentiation is a considerably beneficial model for manufacturers who must meet demands for multiple products that share a common intermediate part. Figure 7 shows the effects of variations of common part's completion rate a on the optimal production cycle time T*. As the completion rate a increases, the optimal cycle time T* decreases significantly, and in the proposed model optimal cycle time T* is reduced by 25.5 % at a = 0.5 (i.e., it decreases from 0.6193 [9] to 0.4614 (years)) as compared to that in prior study which used a single-stage production scheme. Such an analytical result indicates our proposed two-stage multi-item production scheme with delayed differentiation provides a shorter cycle time (or faster response time) than that in a conventional one-stage multi-item system [9]. E\TCU(T, n)\ # 4? / & f / f J J> J J> T E[TCU(T, n)] $2,350,000 $2,300,000 $2,250,000 $2,200,000 $2,050,000 $2,000,000 1-stage production scheme [9] The proposed 2-stage production scheme with linear cost relationship 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95 Fig. 5 The effects of variations of the production cycle time T Fig. 6 The behavior of E[TCU(T, n)] with respect to on the expected system costs E[TCU(T, n)] the common intermediate part's completion rate a — — — 1-stage production scheme [9] _ The proposed 2-stage production scheme with linear cost relationship n ns n 15 n7.>> n 35 n n 55 n B5 n 75 n B5 n 95 Fig. 7 The effects of variations of common part's completion rate a on the optimal production cycle time T* 4.2 Case 2: Analysis of nonlinear relationship of cost relevant variables In this section, we demonstrate that the proposed model is capable of analyzing any given nonlinear relationship between the common part's relevant costs and its completion rate a. For instance, if a nonlinear relationship of 'aA(1/3)' between common part's relevant costs and a is known, then Cq = [aA(1/3)]C1 = [(0.5)a(1/3)]$80 = $63, so it obviously has higher production Advances in Production Engineering & Management 11(4) 2016 341 Chiu, Kuo, Chiu, Hsieh cost (or called value) than that in the linear relationship case (which is $40). Apply the similar computation we have the following values of other relevant parameters: CR,0 = $40, K0 = $13,493, hifl = h4,o = $8, and h2,o = $24. Assume the following parameters' values remain the same as stated in subsection 4.1: Pi,o = 120,000, P20 = 96,000, and xo = [0, 0.04]. Accordingly, in stage 2 we obtain the values of other variables as follows: Cj = $17, $27, $37, $47, and $57, respectively; K = $3,507, $4,007, $4,507, $5,007, and $5,507, respectively; CKi = $10, $15, $20, $25, and $30, respectively; and xi follows the uniform distribution over the intervals [0, 0.01], [0, 0.06], [0, 0.11], [0, 0.16], and [0, 0.21], respectively We apply Eqs. 34, 35, and 25 to obtain the optimal number of shipments n* = 3, the optimal production cycle time T* = 0.4005 (years), and the expected system costs E[TCU(T*, n*)] = $2,093,253. Figure 8 depicts the behavior of E[TCU(T, n)] with respect to the common part's completion rate a under both linear and nonlinear relationships. In nonlinear relationship case, as the common part's completion rate a increases, the expected system costs E[TCU(T, n)] decreases, and it indicates that E[TCU(T, n)] is decreased by 2.45 % at a = 0.5 (i.e., system costs declined from $2,145,834 to $2,093,253) compared to that in the earlier linear case. The analytical results demonstrate that the proposed two-stage multi-item production scheme with delayed differentiation is a greatly beneficial model to manufacturers who have to meet demands for multiple products that share a common intermediate part. Figure 9 illustrates the behavior of the optimal production cycle time T* with respect to the common part's completion rate a under both linear and nonlinear relationships. As completion rate a increases, the optimal production cycle time T* decreases significantly, and in the nonlinear case, the optimal cycle time T* is shortened by 13.20 % at a = 0.5 (i.e., it reduces from 0.4614 to 0.4005) compared to that in the earlier linear case. Therefore, it demonstrates that the proposed two-stage multi-item production scheme with delayed differentiation is a considerably beneficial model (in terms of faster response cycle time) for manufacturers who have to meet demands for multiple products that share a common intermediate part. Furthermore, the analytical results reveals that if the common part's relevant costs are higher (e.g., having a nonlinear relationship aA(1/3) rather than the linear one), then the optimal cycle time T* reduces significantly compared to that in the linear case. E[TCU(T, n)] $2,300,000 $2,250,000 $2,200,000 $2,150,000 $2,100,000 $2,050,000 $2,000,000 1-stage production scheme [9] The proposed 2-stage production scheme with linear cost relationship The proposed 2-stage production scheme with nonlinear cost relationship i i i i i i i i 11 i i i i i i i 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95 — — 1-stage production scheme [9] The proposed 2-stage production scheme with linear cost relationship The proposed 2-stage production scheme with nonlinear cost relationship • —i—i—i—i—i—i—i—i—i—'-i—i—i—i—i—i—i—i—r-0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95 Fig. 8 The behavior of E[TCU(T, n)] with respect to the com- Fig. 9 The behavior of the optimal production cycle mon intermediate part's completion rate a under both linear time T* with respect to the common intermediate and nonlinear relationships part's completion rate a under both linear and nonlinear relationships 5. Conclusion Inspired by the potential benefits derived from applying delayed differentiation to multi-product systems, and with the aim of providing managers of transnational enterprises with information to assist them in achieving the key operational goals such as maximizing machine utilization, ensuring product quality, lowering overall operating costs, and shortening response time, this 342 Advances in Production Engineering & Management 11(4) 2016 Effect of delayed differentiation on a multi-product vendor-buyer integrated inventory system with rework study explores the effect of delayed differentiation on a multi-product vendor-buyer integrated inventory system with rework, using a single machine production scheme. Using mathematical modeling and optimization methods, we derive the closed-form optimal replenishment cycle time and delivery decisions and demonstrate the practical use of our results through a numerical example. The results reveal that our proposed multi-product fabrication scheme with delayed differentiation strategy is considerably beneficial in saving expected system costs and reducing replenishment cycle time. Further analysis also indicates that when the common intermediate part's value is higher (e.g., having a nonlinear relationship aA(1/3) rather than the linear one), both the expected system costs and production cycle time reduces significantly compared to that in the linear case. For future study, to explore and compare the effects of the dual-machine production scheme on the optimal operating policies of the same system would be an interesting direction. Acknowledgement Authors would like to express their gratefulness to the Ministry of Science and Technology of Taiwan for sponsor of this research (under grant no.: MOST 102-2410-H-324-015-MY2). References [1] Hadley, G., Whitin, T.M. (1963). Analysis of Inventory Systems, Prentice-Hall, New Jersey, USA. [2] Silver, E.A., Pyke, D.F., Peterson, R. (1998). Inventory Management and Production Planning and Scheduling, John Wiley & Sons, New York, USA. [3] Nahmias, S. (2009). Production & Operations Analysis, McGraw-Hill, New York, USA. [4] Aggarwal, V. (1984). Grouping multi-item inventory using common cycle periods, European Journal of Operational Research, Vol. 17, No. 3, 369-372, doi: 10.1016/0377-2217(84)90132-2. [5] Rosenblatt, M.J., Rothblum, U.G. (1990). On the single resource capacity problem for multi-item inventory systems, Operations Research, Vol. 38, No. 4, 686-693, doi: 10.1287/opre.38.4.686. [6] Aliyu, M.D.S., Andijani, A.A. (1999). Multi-item-multi-plant inventory control of production systems with short-ages/backorders, International Journal of Systems Science, Vol. 30, No. 5, 533-539, doi: 10.1080/002077299 292272. [7] Balkhi, Z.T., Foul, A. (2009). A multi-item production lot size inventory model with cycle dependent parameters, International Journal of Mathematical Models and Methods in Applied Sciences, Vol. 3, No. 2, 94-104. [8] Rahmani, D., Ramezanian, R., Fattahi, P., Heydari, M. (2013). A robust optimization model for multi-product two-stage capacitated production planning under uncertainty, Applied Mathematical Modelling, Vol. 37, No. 20-21, 8957-8971, doi: 10.1016/j.apm.2013.04.016. [9] Chiu, Y-S.P., Chiang, K-W., Chiu, S.W., Song, M-S. (2016). Simultaneous determination of production and shipment decisions for a multi-product inventory system with a rework process, Advances in Production Engineering & Management, Vol. 11, No. 2, 141-151, doi: 10.14743/apem2016.2.216. [10] Caggiano, K.E., Jackson, P.L., Muckstadt, J.A., Rappold, J.A. (2009). Efficient computation of time-based customer service levels in a multi-item, multi-echelon supply chain: A practical approach for inventory optimization, European Journal of Operational Research, Vol. 199, No. 3, 744-749, doi: 10.1016/j.ejor.2008.08.002. [11] Bjork, K.-M. (2012). A multi-item fuzzy economic production quantity problem with a finite production rate, International Journal of Production Economics, Vol. 135, No. 2, 702-707, doi: 10.1016/i.iipe.2011.10.003. [12] Ma, W.-N., Gong, D.-C., Lin, G.C. (2010). An optimal common production cycle time for imperfect production processes with scrap, Mathematical and Computer Modelling, Vol. 52, No. 5-6, 724-737, doi: 10.1016/j.mcm. 2010.04.024. [13] Guerrero, W.J., Yeung, T.G., Gueret, C. (2013). Joint-optimization of inventory policies on a multi-product multiechelon pharmaceutical system with batching and ordering constraints, European Journal of Operational Research, Vol. 231, No. 1, 98-108, doi: 10.1016/j.ejor.2013.05.030. [14] Chiu, Y.-S.P., Sung, P.-C., Chiu, S.W., Chou, C.-L. (2015). Mathematical modeling of a multi-product EMQ model with an enhanced end items issuing policy and failures in rework, SpringerPlus, Vol. 4, No. 1, 679, doi: 10.1186/ s40064-015-1487-4. [15] Gerchak, Y., Magazine, M.J., Gamble, B.A. (1988). Component commonality with service level requirements, Management Science, Vol. 34, No. 6, 753-760, doi: 10.1287/mnsc.34.6.753. [16] Garg, A., Tang, C.S. (1997). On postponement strategies for product multiple points of differentiation, IIE Transactions, Vol. 29, No. 8, 641-650, doi: 10.1080/07408179708966374. [17] Graman, G.A. (2010). A partial-postponement decision cost model, European Journal of Operational Research, Vol. 201, No. 1, 34-44, doi: 10.1016/j.ejor.2009.03.001. [18] Collier, D.A. (1982). Aggregate safety stock levels and component part commonality, Management Science, Vol. 28, No. 11, 1296-1303, doi: 10.1287/mnsc.28.11.1296. Advances in Production Engineering & Management 11(4) 2016 343 Chiu, Kuo, Chiu, Hsieh [19] Swaminathan, J.M., Tayur, S.R. (1998). Managing broader product lines through delayed differentiation using vanilla boxes, Management Science, Vol. 44, No. 12 (Part 2), S161-S172. [20] Swaminathan, J.M., Tayur, S.R. (2003). Models for supply chains in e-business, Management Science, Vol. 49, No. 10, 1387-1406, doi: 10.1287/mnsc.49.10.1387.17309. [21] Bernstein, F., DeCroix, G.A., Wang, Y. (2011). The impact of demand aggregation through delayed component allocation in an assemble-to-order system, Management Science, Vol. 57, No. 6, 1154-1171, doi: 10.1287/mnsc. 1110.1340. [22] Shih, W. (1980). Optimal inventory policies when stock-outs result from defective products, International Journal of Production Research, Vol. 18, No. 6, 677-686, doi: 10.1080/00207548008919699. [23] Makis, V. (1998). Optimal lot sizing and inspection policy for an EMQ model with imperfect inspections, Naval Research Logistics, Vol. 45, No. 2, 165-186, doi: 10.1002/(SICI)1520-6750(199803)45:2<165::AID-NAV3>3.0.C0: 2-6. [24] Ojha, D., Sarker, B.R., Biswas, P. (2007). An optimal batch size for an imperfect production system with quality assurance and rework, International Journal of Production Research, Vol. 45, No. 14, 3191-3214, doi: 10.1080/ 00207540600711853. [25] Chiu, Y.-S.P., Wu, M.-F., Cheng, F.-T., Hwang, M.-H. (2014). Replenishment lot sizing with failure in rework and an enhanced multi-shipment policy, Journal of Scientific & Industrial Research, Vol. 73, 648-652. [26] Khan, M., Jaber, M.Y., Ahmad, A.-R. (2014). An integrated supply chain model with errors in quality inspection and learning in production, Omega, Vol. 42, No. 1, 16-24, doi: 10.1016/j.omega.2013.02.002. [27] Yeh, T.-M., Pai, F.-Y., Liao, C.-W. (2014) Using a hybrid MCDM methodology to identify critical factors in new product development, Neural Computing and Applications, Vol. 24, No. 3, 957-971, doi: 10.1007/s00521-012-1314-6. [28] Pal, S., Mahapatra, G.S., Samanta, G.P. (2015). A production inventory model for deteriorating item with ramp type demand allowing inflation and shortages under fuzziness, Economic Modelling, Vol. 46, 334-345, doi: 10.1016/j.econmod.2014.12.031. [29] Schwarz, L.B. (1973). A simple continuous review deterministic one-warehouse N-retailer inventory problem, Management Science, Vol. 19, No. 5, 555-566, doi: 10.1287/mnsc.19.5.555. [30] Goyal, S.K., Gupta, Y.P. (1989). Integrated inventory models: The buyer-vendor coordination. European Journal of Operational Research, Vol. 41, No. 3, 261-269, doi: 10.1016/0377-2217(89)90247-6. [31] Hill, R.M. (1996). Optimizing a production system with a fixed delivery schedule, Journal of the Operational Research Society, Vol. 47, No. 7, 954-960, doi: 10.1057/jors.1996.121. [32] Thomas, D.J., Hackman, S.T. (2003). A committed delivery strategy with fixed frequency and quantity, European Journal of Operational Research, Vol. 148, No. 2, 363-373, doi: 10.1016/S0377-2217(02)00398-3. [33] £omez, n., Stecke, K.E., £akanyildirim, M. (2012). Multiple in-cycle transshipments with positive delivery times, Production and Operations Management, Vol. 21, No. 2, 378-395, doi: 10.1111/j.1937-5956.2011.01244.x. [34] Wu, M.-F., Chiu, Y.-S.P., Sung, P.-C. (2014). Optimization of a multi-product EPQ model with scrap and an improved multi-delivery policy, Journal of Engineering Research, Vol. 2, No. 4, 103-118, doi: 10.7603/s40632-014-0027-7. [35] Safaei, M. (2014). An integrated multi-objective model for allocating the limited sources in a multiple multi-stage lean supply chain, Economic Modelling, Vol. 37, 224-237, doi: 10.1016/j.econmod.2013.10.018. [36] Chiu, S.W., Sung, P.-C., Tseng, C.-T., Chiu, Y.-S.P. (2015). Multi-product FPR model with rework and multi-shipment policy resolved by algebraic approach,Journal of Scientific & Industrial Research, Vol. 74, No. 10, 555-559. [37] Tseng, C.-T., Wu, M.-F., Lin, H.-D., Chiu,Y.-S.P. (2014). Solving a vendor-buyer integrated problem with rework and a specific multi-delivery policy by a two-phase algebraic approach, Economic Modelling, Vol. 36, 30-36, doi: 10.1016/j.econmod.2013.09.013. [38] Hishamuddin, H., Sarker, R.A., Essam, D. (2014). A recovery mechanism for a two echelon supply chain system under supply disruption, Economic Modelling, Vol. 38, 555-563, doi: 10.1016/j.econmod.2014.02.004. [39] Chiu, S.W., Huang, C.-C., Chiang, K.-W., Wu, M.-F. (2015). On intra-supply chain system with an improved distribution plan, multiple sales locations and quality assurance. SpringerPlus, Vol. 4, No. 1, 687, doi: 10.1186/ s40064-015-1498-1. [40] Pai, F.-Y. (2015). How supplier exercised power affects the cooperative climate, trust and commitment in buyer-supplier relationships, International Journal of Business Excellence, Vol. 8, No. 5, 662-673, doi: 10.1504 /IIBEX.2015.071275. [41] Ocampo, L.A. (2015). A hierarchical framework for index computation in sustainable manufacturing, Advances in Production Engineering & Management, Vol. 10, No. 1, 40-50, doi: 10.14743/apem2015.1.191. [42] Saric, T., Simunovic, K., Pezer, D., Simunovic, G. (2014). Inventory classification using multi-criteria ABC analysis, neural networks and cluster analysis, Technical Gazette - Tehnicki vjesnik, Vol. 21, No. 5, 1109-1115. [43] Kaveh, M., Dalfard, V.M. (2012). A study on the effect of inflation and time value of money on lot sizing in spite of reworking in an inventory control model, Technical Gazette - Tehnicki vjesnik, Vol. 19, No. 4, 819-826. [44] Al-Hawari, T., Ahmed, A., Khrais, S., Mumani, A. (2013). Impact of assignment, inventory policies and demand patterns on supply chain performance, International Journal of Simulation Modelling, Vol. 12, No. 3, 164-177, doi: 10.2507/HSIMM12(3)3.235 [45] Rardin, R.L. (1998). Optimization in Operations Research, Prentice-Hall, New Jersey, USA. 344 Advances in Production Engineering & Management 11(4) 2016 APEM jowatal Advances in Production Engineering & Management Volume 11 | Number 4 | December 2016 | pp 345-354 http://dx.doi.Org/10.14743/apem2016.4.232 ISSN 1854-6250 Journal home: apem-journal.org Original scientific paper Experimental modeling of fluid pressure during hydroforming of welded plates Karabegovic, E.a*, Poljak, J.a aFaculty of Technical Engineering, University of Bihac, Bihac, Bosnia and Herzegovina A B S T R A C T A R T I C L E I N F O The procedure of hydroforming belongs to one of the modern methods of sheet and tube design, usually of complex configuration. Research in the field of plastic forming using fluids usually relates to the analysis of important parameters that would enable high-quality design of elements and execution of the process in stable conditions. The hydroforming process of welded sheets found its application in manufacturing of tanks and other sheet parts in automotive industry, where, in addition to technical and technological characteristics of the obtained piece, it is necessary to achieve stability of the process and its economic feasibility. Experimental research in this paper had been aimed at the analysis of results and modeling of working fluid pressure during hydroforming of welded sheets of two kinds of material (St 37 and Al 99.5) for two sheet thicknesses (1.5 mm and 2.0 mm). Modeling was done by regression method, whose analysis is the determination of functional relationships between a dependent variable and two independent variables. Application of mathematical modeling method enabled working fluid pressure which confirmed the impact of input variables of hydroforming process (yield strength and sheets thickness) onto the values of working fluid pressure. Experimental results obtained for working fluid pressure enabled easier planning and projection of hydroforming process. © 2016 PEI, University of Maribor. All rights reserved. Keywords: Forming Hydroforming Welding sheet metal Fluid pressure Modelling Regression *Corresponding author: edina-karabeg@hotmail.com (Karabegovic, E.) Article history: Received 15 September 2016 Revised 20 November 2016 Accepted 24 November 2016 1. Introduction The development of the automotive industry is based on the emergence of new materials and technologies for their processing. New methods of processing materials enable the achievement of technical and technological characteristics required by the market of the finished product. Hydroforming is the process of forming sheets and tubes during which, by function of fluids underpressure, pieces, most often of complex shapes, are formed for the automotive industry [1-4]. The process is quick, inexpensive and meets the quality of shaped elements. So far, a number of studies have been conducted that included analytical, experimental, numerical, mathematical and other analyses of various treatment processes [5-10]. Hein and Vollertsen (1999) have conducted experimental and numerical research in order to establish the technological and economic characteristics of the process during hydroforming of double sheets [11]. The analysis and comparison of conventional deep drawing and process of sheet element hydroforming by using the finite element method [12] is the area of research of Chang et al. (2001). By application of CCD (Charge Coupled Device) camera, Groche, P. et al. (2007) conducted the control of sealing system during the execution of sheet hydroforming process. 345 Karabegovic, Poljak The control of fluid pressure upon hydroforming of double sheets is an area of research of Assempour and Emami (2009) [13]. In their paper Ertugrul et al. (2009) analyzed hydroforming of laser welded sheets [14]. In addition, Liewald and Bolay (2010) in their paper state the analysis of hydroforming process of double sheets [15]. By application of finite element method (FEM) Zhang et al. (2015) analyzed the impact of stress or pressure of fluids on the improvement of quality of the double sheets hydroforming [16]. As can be seen, there is a constant need to improve the technical conditions and methods for execution of the hydroforming process. This paper gives the analysis of the experimental execution of hydroforming process of welded sheets of two types of materials with different thicknesses. During hydroforming of pieces with defined size and shape, we have measured the working fluid pressure and displacement (expansion) of welded sheets. The experimental results were used in the process of mathematical modeling. Mathematical model for working fluid pressure was obtained for defined conditions of process execution. 2. The hydroforming process Hydroforming of connected sheets is a technique of forming by effect of fluid pressure in the interior of the welded sheets. Thereby, the sheets are deformed (spread) to the interior shape of the die which defines the shape and dimensions of the finished piece. The hydroforming process of welded sheets is used for making fuel tanks, the car doors, etc. Execution process scheme is given in four phases, as shown in Fig. 1: a) Placement of welded sheets into a die (matrix), b) Activation of pressure and pre-forming, c) Calibration with final pressure of the fluid, d) The finished piece is removed from the die. With this process it is possible to form welded sheets of the same or different material and thickness. For the analysis in this paper, we conducted hydroforming of work element of the defined shape. (a) (b) final of forming m i i finished piece M fd) Fig. 1 Scheme of welded sheet hydroforming process 3. Experimental measurement of working fluid pressure The following parameters were defined for the execution of the experiment: • Geometrical shape and dimensions of finished pieces, • Material (sheet) to produce beginning piece, • Method of sheet welding, • Tools for hydroforming and fluid, 346 Advances in Production Engineering & Management 11(4) 2016 Experimental modeling of fluid pressure during hydroforming of welded plates • Pressure system or measuring amplifier device (pump), • System for control and measurement of process parameters: computer and measurement equipment for information gathering. Fig. 2 gives geometric shape and dimensions of the finished piece to be formed [17]. v/ Fig. 2 Geometric shape and dimension of finished piece 3.1 Beginning piece and tools for hydroforming When producing pieces (beginning piece), two types of materials (sheets) were applied, aluminum (99.5 % ) and steel (St 37). Raw parts were produced from cold rolled sheets with thickness: s = 1.5 mm and s = 2.0 mm. Mechanical properties of the material for raw parts are given in Table 1. Table 1 Mechanical properties of materials Types of materials Yield strength, N/mm2 Mechanical strength, N/mm2 Modulus of elasticity, N/mm2 Steel St 37 235 410 2.1106 Aluminium 99.5 100 120 0.7 ■ 104 MIG welding procedure (protective argon gas) was selected for welding of the two materials. Fig. 3 depicts the position of raw part in the tool for hydroforming of welded sheets. Tool for hydroforming of welded sheets was produced from structural steel St 37 and consists of two parts (upper and lower matrix) and connection bolts. Fig. 3 The position of beginning piece in the tool for hydroforming of welded sheets 3.2. Pump for achieving working fluid pressure Hydraulic high-pressure pump pmax = 3-107 Pa was used for the execution of the experiment Working fluid for hydroforming is oil „Inol hidrol-X46", of density 0.884 g/cm3 (20 °C). Advances in Production Engineering & Management 11(4) 2016 347 Karabegovic, Poljak 3.3 Measuring equipment Measuring amplifier device "Spider 8" with eight independent measuring channels was used for the measurement of working fluid pressure and displacement (expansion) of sheets, as shown in Fig. 4. measuring amplifier device "Spider 8" pressure sensor Fig. 4 Measuring equipment sensor feed Applied sensor for the measurement of the working pressure-P8AP, given in Fig. 5, operates on the principle of strain gauges. The nominal sensitivity of the sensor is 2 mV/5-107 Pa. Sensor measuring range is 0-5-107 Pa. Displacement sensor-WA20, shown in Fig. 4, operates on inductive principle. Nominal sensitivity of the sensor is 80 mV/20 mm. Measuring range of the sensor is 0-20 mm. Fig. 5 depicts the position of the sensor for measuring fluid pressure and displacement (expansion) of sheets during hydroforming of welded sheets. Î sensor feed pressure sensor beginning piece ' H CI Fig. 5 The position of sensors during hydroforming of welded sheets and shaped piece 3.4 Number of tests in the experiment Experimental measurements of the defined parameters are aimed at practical application of the results obtained in the planning, design and implementation process. Research and analysis of processes are important for achieving greater process stability in the given circumstances. Number of probes in this experiment was determined by the rotatable plan of the experiment and expression: N= 2li~p + 2k + n0 =nk +na +n0 = 22 + 2-2 + 5 = 13 (1) tool 48 Advances in Production Engineering & Management 11(4) 2016 Experimental modeling of fluid pressure during hydroforming of welded plates where N is total number of experiments, nk is number of change of variables, no is number of repetitions in plan center, and na is number of symmetrically positioned points at plan center. Experiment plan matrix has the shape presented in Table 2. Experiment was conducted in the laboratory of the Faculty of Technical Engineering of the University of Bihac. Table 2 Experiment plan matrix Input variables of the process Physical values Coded values Trial No. Yield strength, N/mm2 (0D.2) Thickness, mm ( si) Xo Xi X2 Xi X2 Xi2 X22 Vector output Yi 1 100 1.5 1 -1 -1 1 1 1 Yi 2 235 1.5 1 1 -1 -1 1 1 Y2 3 100 2.0 1 -1 1 -1 1 1 Y3 4 235 2.0 1 1 1 1 1 1 Y4 5 168 1.75 1 0 0 0 0 0 Y5 6 168 1.75 1 0 0 0 0 0 Y6 7 168 1.75 1 0 0 0 0 0 Y7 8 168 1.75 1 0 0 0 0 0 Ys 9 168 1.75 1 0 0 0 0 0 Y9 10 263 1.75 1 1.414 0 0 2.0 0 Y10 11 72 1.75 1 -1.414 0 0 2.0 0 Y11 12 168 2.10 1 0 1.414 0 0 2.0 Y12 13 168 1.39 1 0 -1.414 0 0 2.0 Y13 4. Measurement results Experimental results of working fluid pressure and displacement (expansion) of sheets during hydroforming of welded aluminum and steel sheets, with thickness of 1.5 mm and 2.0 mm, are given in Fig. 6. Experimental results of working fluid pressure for hydroforming of welded aluminum sheets are given in Fig. 7. Experimental results of working fluid pressure for hydroforming of welded steel St 37 sheets are given in Fig. 8. Comparative results of the working fluid pressure for hydroforming of welded sheets of aluminum and steel with 1.5 mm thickness are given in Fig. 9. Comparative results of the working fluid pressure for hydroforming of welded sheets of aluminum and steel with 2.0 mm thickness are given in Fig. 10. Fig. 6 Experimental results of working fluid pressure and displacement during hydroforming of welded sheets Advances in Production Engineering & Management 11(4) 2016 349 Karabegovic, Poljak Working fluid pressure for hydroforming of welded Al sheets Number of probe Fig. 7 Working fluid pressure for hydroforming of welded aluminum sheets Working fluid pressure for hydroforming of welded St 37 sheets CN ^ 100 2 « 80 o tf 60 S- i*-000000^ Steel 5=1.5 mm Steel 5=2.0 mm £ 40 i-, a. 12 3 Number of probe Fig. 8 Working fluid pressure for hydroforming of welded steel sheets Working fluid pressure for hydroforming Al and St 37 5=1.5 mm Number of probe Fig. 9 Comparative results of working fluid pressure for aluminum and steel sheets with 1.5 mm thickness N 100 B 2 80 S 60 Č 40 m P 20 Working fluid pressure for hydroforming Al and St 37 5= 2.0 mm -Al -Steel Number of probe 350 Fig. 10 Comparative results of working fluid pressure for aluminum and steel with 2.0 mm thickness Advances in Production Engineering & Management 11(4) 2016 Experimental modeling of fluid pressure during hydroforming of welded plates Analysis of the experiment results The analysis of the obtained experimental values of the working fluid pressure for 13 samples (welded sheet 1.5 mm and 2.0 mm thickness, aluminum Al 99.5 % and steel St 37) provides the following conclusions: • Average value of the working fluid pressure to form 2.0 mm thick aluminum sheet compared to the sheet thickness of 1.5 mm is higher for around 10.38 %, • Average value of the working fluid pressure to form 2.0 mm thick steel sheet compared to the sheet thickness of 1.5 mm is higher for around 7.67 %, • Increase in the working fluid pressure in both materials is the result of an increase in the thickness of welded sheets, • The design of welded 1.5 mm thick steel sheet requires higher values of the working fluid pressure for around 42.84 % compared to welded aluminum sheets of the same thickness, • The design of welded 2.0 mm thick steel sheet requires higher values of the working fluid pressure for around 41.12 % compared to welded aluminum sheets of the same thickness, • Increase in the working fluid pressure during hydroforming of welded steel sheets is justified due to the differences in the mechanical properties of steel and aluminum, • Average value of displacement deviations for 13 samples is about 7 %, which is caused by changes in the structure and quality of the sheet. 5. Modeling of working fluid pressure during hydroforming od welded sheets The measured values of the working fluid pressure and displacement after conducted 13 experiments will be used to define a mathematical model that, in the appropriate level of accuracy, adequately describes the hydroforming process of welded sheets with defined shape and conditions of the process execution. The model would be used in the design phase of the process, in analyzing and forecasting the state of the process [18]. The analysis of the stochastic process starts from the general functional relationship between the dependent variable size (7,) and independent variables (x,), which can be presented with the model as follows: Yi f(Xi) f{x1,x2) i.e. Y = p = f(ao.2,s) Coded values of physical quantities are obtained using the following expressions: (2) (3) = 1 + 2 ln< 0.2; In a 0.2; lnah ? —ln5 then Jidoka or quality on source, should be applied. With all three indexes defined, RPN can be calculated as the end of the plan phase. After the plan phase, do phase or realization phase follow. For all defined RPN, corrective actions or Kaizen, should be taken. Suggested improvements had to be set on a list of solutions, with exact deadlines and with responsible team members. Third phase is check phase where the action taken has to be checked with repeating of RPN calculation. If there is no progress, then Kaizen should be performed again. This check phase have to be realized very carefully, because after this phase failure should be standardized. Last phase is act phase. When the solution is finally found and the progress has been made, failure and elected solution have to be standardized and ready to use. 5. The new PFMEA approach: Case study from automotive industry The company elected for the case study is automotive company which produce electronic circuit boards and electronic cables for automobiles. Company is supplying automotive suppliers and corporations all over the world. This company applies lean approach in their production system for a long time and also use PFMEA for prevention of failures and risks. The use of PFMEA in this company is obligatory according to ISO9000 and ISO TS16949 standards. PFMEA for product MSM6BL was already realized on traditional way. Results are presented in Table 2. Measured conditions taken for comparation are: number of team members, identification of failures, change of S, O, D and RPN values, reduced S, O, D, and RPN values after taken actions, RPN value over 100, and S, O, D value over 8. These conditions are measured in total amount for whole PFMEA. The goal is to compare them with the new approach and see the differences after its implementation. S, O, D and RPN value changes are also calculated in total change regardless if it is increased or reduced value. Methodology set in algorithm form from Fig.1 is used for realization of PFMEA. Realized PFMEA report is shown in Appendix 1. The changes made after the new approach are painted in grey in Appendix 1. They have been implemented and used for calculation. The data for a new state are presented in Table 2. These data were taken from Appendix 1, also. After this step, the comparation between the state after traditional approach to PFMEA realization and the new approach to PFMEA realization with integrated lean approach, was made. These results are also presented in Table 2. From Table 2, it can be seen that almost all conditions are changed, except one - S. Two conditions; the actions taken and identification of failures, are very important for analysis because they are related to the purpose of analysis to detect a failure and take action to improve it. Gen- 360 Advances in Production Engineering & Management 11(4) 2016 An integrated lean approach to Process Failure Mode and Effect Analysis (PFMEA): A case study from automotive industry chi Genbutsu have mostly contributed to these changes. Lean approach stimulated involvement of more employees, including six workers. As it was predicted, workers contributed a lot to failure identification because they are directly involved in a process. But Genchi Genbutsu stimulated PFMEA team to go directly in process and observe an actual status. Increased identification of failures and the taken action can avoid hidden failures reaching customer. Also, they affected the changes of S, O, D indexes and RPN values. Moreover, during these changes, value of some S, O, D indexes exceeded 8 and RPN values exceeded 100. The lack of situations when S, O, D and RPN exceed predicted values for improvements, may cause problems if failures reach a customer. Only one condition which is not reduced is S index due to the actions taken, which is not a big issue due to the O and D indexes reduction for that failures. This was achieved due to the application of Kaizen. Also, some of D indexes were reduced due to the application of Jidoka. Causes of failures were superficially defined and some resulted with mixing of causes and effects. With application of lean tool - 5 why, root causes of failures were deeply analysed. The actions taken were oriented on fixing root causes of failure, not effects. One more lean tool used in this case study, is standardized work. Some of failures were standardized. This means that in the next PFMEA for some of the new processes or products, standardized solution will be used. That will save a lot of time. Along these integrated lean tools into PFMEA, Poka-Yoke and 5S were also used as the lean tool for recommended actions during PFMEA realization. Table 2 Comparison of state before and after lean approach integration in PFMEA Measured conditions State after traditional approach New state after lean approach Improvements (%) Number of team members 2 9+6 workers 85 The actions taken 1 19 95 Identification of failures 18 27 33 Change of S value 155 227 32 Change of O value 64 105 39 Change of D value 85 161 47 Change of RPN value 1642 3720 56 Reduced RPN value due to taken actions 2602 3366 23 Reduced S value due to taken actions 226 226 / Reduced O value due to taken actions 80 127 37 Reduced D value due to taken actions 121 139 13 RPN value over 100 / 16 100 S, O, D values over 8 / 3 100 6. Conclusion The new approach with the integration of lean approach into PFMEA for improvement of specific shortcomings, is presented in this paper. This new PFMEA approach has proven to be very good and practical combination for identifying and fixing problems and failures. The case study realized from the automotive industry was used for the new approach testing. The state with traditional PFMEA is compared with the new state where almost all measured conditions were changed and improved. Therefore, this approach was practically approved. The main advantage of this approach was in improvement of process and product quality, which is mostly important for customers. For a change, this approach is applicable in real-world cases in every process or industry where the lean approach is implemented. Also, it is very simple and practical to use and does not require big investments, implementation of new technology and complicated additional education. The way this approach has stopped potential failures to reach a customer was not identified with traditional approach. This was presented on practical example. Therefore, this new approach increased identification and prevention possibilities as well. Advances in Production Engineering & Management 11(4) 2016 361 Banduka, Veza, Bilic This new approach have also few constraints. It is used for specific group of shortcomings and it is not comprehensive. Also cannot be applicable in processes and industries where the lean approach has not been implemented. Costs are not included and there is different aspect of looking on it from perspective of PFMEA and lean approach. Future work should be oriented to fixing all other shortcomings of PFMEA, costs especially. Lean approach is oriented to long term thinking about costs. From PFMEA team perspective is different. During PFMEA realization team have to percept immediately if the solution is profitable. Thus in further research, the balance should be found between urgent need for costs from PFMEA aspect and long term thinking profitability from the lean approach aspect. Acknowledgement This paper has been supported by Croatian Science Foundation under the project Innovative Smart Enterprise -INSENT (1353). This publication also has been supported by the European Commission under the Erasmus Mundus project Green-Tech-WB:Smart and Green technologies for innovative and sustainable societies in Western Balkans (551984-EM-1-2014-1-ES-ERA MUNDUS-EMA2). References [1] Harrison, A., van Hoek, R. (2005). Logistics management and strategy, Pearson Education, Harlow, England. [2] Stamatis, D.H. (2003). Failure mode and effect analysis: FMEA from theory to execution, ASQ Quality Press, Milwaukee, Wisconsin, USA. [3] ISO 9000. Quality management systems - Fundamentals and vocabulary, FMEA obligations, from http://www.bti.secna.ru/education/smq/docs/news/iso 9000 2005 e.pdf, accessed October 30, 2016. [4] ISO/TS 16949:2002(E). Quality management systems - Particular requirements for the application of ISO 9001: 2000 for automotive production and relevant service part organizations, FMEA obligations, from https://ciiaas.files.wordpress.com/2007/11/iso-ts-16949-2002.pdf. accessed June 13, 2016. [5] Johnson, K.G., Khan, M.K. (2003). A study into the use of the process failure mode and effects analysis (PFMEA) in the automotive industry in the UK, Journal of Materials Processing Technology, Vol. 139, No. 1-3, 348-356, doi: 10.1016/S0924-0136(03)00542-9. [6] Liu, H.C., Liu, L., Liu, N. (2013). Risk evaluation approaches in failure mode and effects analysis: A literature review, Expert systems with applications, Vol. 40, No. 2, 828-838, doi: 10.1016/j.eswa.2012.08.010. [7] Biddle, J. The lean benchmark report, from http://consumergoods.edgl.com/column/The-Lean-Benchmark-Report48987. accessed June 13, 2016. [8] Lorenc, M., Jentzsch, A., Andersen, M., Noack, B., Waffenschmidt, L., Schuh, G., Rudolf, S. The lean advantage in engineering: developing better products faster and more efficiently, from https://www.bcgperspectives.com, accessed June 13, 2016. [9] Chrysler LLC, Ford Motor Company, General Motors Corporation, Potential failure mode and effect analysis (FMEA), Reference Manual, 4th edition, from www.engmatl.com/home/finish/20-engineering.../160-fmea-manual. accessed June 13, 2016. [10] Korenko, M., Krocko, V., Kaplík, P. (2012). Use of FMEA method in manufacturing organization, Journal of Manufacturing and Industrial Engineering, Vol. 11 No. 2, 48-50. [11] Montgomery, T.A., Pugh, D.R., Leedham, S.T., Twitchett, S.R. (1996). FMEA automation for the complete design process, In: Proceedings of Reliability and Maintainability Symposium - International Symposium on Product Quality and Integrity, IEEE, Las Vegas, NV, USA, 30-36. [12] Vujica-Herzog, N., Tonchia, S. (2014). An instrument for measuring the degree of lean implementation in manufacturing, Strojniški vestnik - Journal of Mechanical Engineering, Vol. 60, No. 12, 797-803, doi: 10.5545/sv-jme.2014.1873. [13] Liker, J.K. (2004). The Toyota way, 14 management principles from the world's freatest manufacturer, McGraw-Hill, New York, USA. [14] Womack, J.P., Jones, D.T. (2010). Lean thinking: Banish waste and create wealth in your corporation, Simon and Schuster, New York, USA. [15] Shekari, A., Fallahian, S. (2007). Improvement of lean methodology with FMEA, In: Proceedings of the 18th Annual Conference - POMS (Production and Operation Management Society), Texas, USA, 007-0520. [16] Sawhney, R., Subburaman, K., Sonntag, C., Rao Venkateswara Rao, P., Capizzi, C. (2010). A modified FMEA approach to enhance reliability of lean systems, International Journal of Quality & Reliability Management, Vol. 27, No. 7, 832-855, doi: 10.1108/02656711011062417. [17] Shahrabi, M., Shojaei, A.A. (2014). Application of FMEA and AHP in lean maintenance, International journal of Modern Engineering Science, Vol. 3, No. 1, 61-73. [18] Puvanasvaran, A.P., Jamibollah, N., Norazlin, N., Adibah, R. (2014). Poka-Yoke Integration into process FMEA, Australian Journal of Basic and Applied Sciences, Vol. 8, No. 7, 66-73. 362 Advances in Production Engineering & Management 11(4) 2016 An integrated lean approach to Process Failure Mode and Effect Analysis (PFMEA): A case study from automotive industry [19] Puvanasvaran, A.P., Jamibollah, N., Norazlin, N. (2014). Integration of POKA YOKE into process failure mode and effect analysis: A case study, American Journal of Applied Sciences, Vol. 11 No. 8, 1332-1342, doi: 10.3844/ ajassp.2014.1332.1342. [20] Teoh, P.C., Case, K. (2007). An evaluation of failure modes and effects analysis generation method for conceptual design, International Journal of Computer Integrated Manufacturing, Vol. 18, No. 4, 279-293, doi: 10.1080/ 0951192042000273122. [21] Arabian-Hoseynabadi, H., Oraee, H., Tavner, P.J. (2010). Failure modes and effects analysis (FMEA) for wind turbines, International Journal of Electrical Power & Energy Systems, Vol. 32 No. 7, 817-824, doi: 10.1016/j.iiepes. 2010.01.019. [22] Moore, C. Failure Modes and Effects Analyses (FMEA) and critical items list (CIL) for the ECS project; EOSDIS core system project, from https://prod.nais.nasa.gov/eps/eps data/161793-AMEND-002-002.pdf. accessed November 13, 2016. [23] Liker, J., Convis, G.L. (2011). The Toyota way to lean leadership: Achieving and sustaining excellence through leadership development, McGraw Hill Professional, New York, USA. Appendix 1 FMEA Number ] Part Name/Part No| Realized by: 51 BSH MSM6BL E index: A1 Project manager FMEA Date (Orig) | FMEA Revision Date Prepared by: 2/18/2015 20.04,2015. Rew Project Leader FMEA Team Project Leader x2, Pro er. Logistics manager, Procurement manager, R&Dx2, Maintananee manager + 6 workers 01 Q C Action Results Process c E Ol Potential Poten ¡al Effect (s) of Failure > o u Potential Cause(s) of Failure c Ol Current Process Current Process Detection z Recomended Responsibility 5 No./Function o- dl Mode i/i IV □ O Prevention Controls: Detaction EC Action(s) E O u č >- c o Ol ■ a H Actions Taken « =c 3 O ? 'S a a. E110 Axial inser ion E110 Axial Component PCB with No Intern Process No insertion failed failed components schooling, approved 72 function; No in storage * working function instruction E110 Axial insertion Component with false value Badly equipped PCB; PCB with failed 6 Wrong components in storage 2 Intern schooling, working instruction Process approved 6 72 No 0 E110 Axial insertion Wrong orientation of Incorect or no function 6 Machine failure 2 Intern schooling, working Process approved 6 72 No 0 E412 Manual inserting E412 Manual Mechanical Trouble in Inattention Intern Visual control; Increase Project Work inserting damaged process 5 schooling, Process 120 workers manager instructions 30 component working instruction approved precaution E412 Manual Component PCB with Inattention Visual Visual Poka yoke Project Aditional tool inserting failed on failed inspection inspection by 125 manager for chack out 30 PCB function, no function before use end inspection of components E412 Manual inserting Wrong component PCB with false components 7 Inattention 4 Intern schooling, working instruction Process approved 8 224 kaizen workshop Project manager Process improvement + check out after component 7 1 8 56 PFMEA report 1/4 Advances in Production Engineering & Management 11(4) 2016 363 Banduka, Veza, Bilic Process No./Function Requirements | Potential Failure Mode Potenial Effect(s) of Failure Severity | Classification | Potential Cause(s) of Failure Occurrence | Current Process Controls: Prevention Current Process Controls: Detaction Detection | RPN Recomended Action(s) Responsibility Target Completetion Date | Action Results Actions Taken tnective i natp 1 1 Occurrence 1 1 Severity 1 Detection | RPN E412 Manual inserting Contact badly placed Trouble in process 8 Inattention 4 Intern schooling, working instruction Visual control; Process approved 8 256 Kaizen + poka yoke Project manager Process improvement + tool for contact set up 8 1 8 64 E500 Wave soldering E500 Wave Soldering Too much solder on solder joint Short circuit onPCB 7 False setting the process parameter -speed of soldering line is wrnnp 4 Intern schooling, working instruction Visual inspection, check the parameter before use 5 140 Preventive machine maintanance Maintanance manager Check out of parameters before start of machine 7 1 4 28 E500 Wave Soldering Use of false temperatur e profil Bad or not soldered joints between pads and component 8 § Inattention, false setting of machine parameters 5 Intern schooling, working instruction Visual inspection 5 200 Preventive machine maintanance Maintanance manager Check out of parameters before start of machine 8 1 6 48 E500 Wave Soldering Uncorect solder joint Unreliable joint -unreliable function 8 False setting the process parameter - 5 Intern schooling, working instruction Visual inspection, check the parameter before use 5 200 Preventive machine maintanance Maintanance manager Check out of parameters before start of machine 8 1 6 48 E500 Wave Soldering Not enough solder on solder joint Bad solder joint 8 False setting the process parameter 4 Intern schooling, working instruction Visual inspection -check the parameter 5 160 Preventive machine maintanance Maintanance manager Check out of parameters before start of machine 8 1 5 40 E500 Wave Soldering Cold contact Bad solder joint 7 False setting the process parameter 4 Intern schooling, working instruction Visual inspection, check the parameter before use 5 140 Preventive machine maintanance Maintanance manager Check out of parameters before start of machine 7 1 5 35 E650 Function test E650 Function test Skipped the final control Bad parts can be sent to the customer 8 Inattention 3 Intern schooling, working instruction Process approved 8 192 One piece flow Worker Work instructions 7 1 5 35 E521 Rework; visual inspection Zabranjeno je neovlaščeno kopiranje i prenošenje prava trečim licima! PFMEA report 2/4 Process No./Function Requirements | Potential Failure Mode Potenial Effect(s) of Failure Severity | Classification | Potential Cause(s) of Failure Occurrence | Current Process Controls: Prevention Current Process Controls: Detaction Detection | RPN Recomended Action(s) Responsibility Target Completetion Date | Action Results Actions Taken (4) where: Advances in Production Engineering & Management 11(4) 2016 369 Klancnik, Hrelja, Balic, Brezocnik a"(t) - Gravitational acceleration of body a in the n-dimension space ^vTCO - Force with stochastic characteristic in specific time t Mia(t) - Inertia mass of body a in specific time t On the basis of the updated equation for the force with stochastic characteristic in specific time t, and on the basis of the gravitational acceleration of body a in the d-dimension space, the new velocity of particle or body in the space with n-dimensions - Eq. 5, can be expressed; the same also applies to the newly calculated position of the body in the solution space - Eq. 6: v2(t + 1) = randa^v2(t) + aW) (5) z£(t + 1)=z£(t) + i£(t + 1) (6) Terms in Eq. 5 and 6: + 1) - New velocity of body a in newly calculated time t + 1 x%(t + 1) - New position of body a in newly calculated time t + 1 randa - Random number at interval [0,1] ensures randomization va (0 - Velocity of motion of body a in specific time t (t) - Acceleration of body a in specific time t xa (0 - Position of body a in specific time t n - Index n stands for the n-dimension space with solutions Gravitational and inertia masses of bodies in the system can also be written as a fitness function to define the success of the algorithm in searching for the solution. It must be borne in mind that larger masses are better indicators of solutions and are more effective; however, as those bodies have larger mass, they move slower in the space with solutions. That phenomenon is compensated for by a random number of lighter bodies representing worse solutions but, in combination with heavier bodies, they take part in converging to proper solutions much faster. If it is assumed that gravitational and inertia masses are identical, they can be written in the form of an equation. Updated equations of the mass system are: Man = Mvn =Min = Mn.n = 1,2.....N (7) = fita(t) -worstjt) ma{t) best(t) —worst(t) (8) (9) where fita(t) stands for the fitness value of body a in specific time t, while best(t) and worst(t) are values representing the best and worst current values of results and can be written even in more detail, particularly, when searching for minimums and maximums. For minimum: best(t) = minje{li2.....N)fitb(t) worst(t) = maXje(li2i_iN)fitb(t') ( ) For maximum: bestit) =maxJe(i.....N)fitb(t) worst(t) = min]e(1,,Nyfitb(t) ( ) To avoid efficiently the local minimums at the beginning of the algorithm start-up the research principle must be introduced, whereas, after n-iterations, the research limit is dimmed by means of elimination and the system of algorithm functioning passes into the search mode. That is reached by gradual elimination and each subsequent algorithm iteration assures better convergence of solutions. Each subsequent potential optimal solution can be designated N-best, where each iteration changes the action of the virtual gravitational force with which it affects all the remaining bodies in the space of solutions. Consequently, N-best is the time function having 370 Advances in Production Engineering & Management 11(4) 2016 Multi-objective optimization of the turning process using a Gravitational Search Algorithm (GSA) and NSGA-II approach the value No at the beginning of time t. In that way, at the beginning of the algorithm, all elements act with gravitational force on all remaining bodies, while the force decreases linearly with the reduction of the number of objects in the system and/or convergence of solutions to one point, which can also be described by means of gravitational collision of two bodies. The algorithm finalises the optimization, when only one body is still available with the greatest possible mass with which it acted on the remaining bodies and attracted them gravitationally. In that way, the Eq. 3 can be transformed and the result is: TO) = ^ randbF2b(t) (12) bEN-best.b^a Fa(.t) - Force with stochastic characteristic in specific time t randb - Random number at interval [0,1], it assures randomization Fab (0 - Force which is caused by the b-th mass particle and acts on the a-th mass particle in time t N-best - Set of first bodies with largest mass and highest fitness value The pseudocode of the GSA algorithm is presented in Fig. 1. 1: Start GSA 2: Set-up search space 3: For each particle 4: Initialize mass particle 5: Start regression analysis module/with measured experimental data 6: END 7: Do 8: For each particle/import measured experimental data 9: Calculate capacity/fitness function of particle 10: Capacity of particle > best capacity of particle (best (t)) ^ new best 11: Evaluate worst values as worst (t) 12: END 13: For each particle 14: Calculate velocity, masses and accelerations of mass particle 15: Update particle position and velocity 16: END 17: While ^ number N of iterations is reached 18: END GSA Fig. 1 Pseudocode of the GSA algorithm 2.3 NSGA-II algorithm background The NSGA-II (elitist Non-dominated Sorting Genetic Algorithm) is a genetic algorithm designed for multi-objective optimization [26]. The algorithm uses the crowding distance metric and non-dominated sorting as its main features. Fig. 2 presents the algorithm's pseudocode. Roughly speaking, an algorithm consists of initialization, selection and re-combination. Selection and re-combination are the main algorithm steps (see Fig. 2). In the first step, the old population of parents and descendants is combined into a joint population of size 2n. In the next step, the subjects from that combined population Rt-1 are sorted front by front using non-dominated sorting. Subjects not dominated by any other subject go into the first front. The first i fronts still going whole into the new parent population Pt are written into that population. The first next front (front i+1) not going whole into the newly created population is sorted by using the crowding metric. The least crowded subjects from front i+1 are added to the new population Pt. The solutions have the largest range if the extreme subjects in the front are evaluated as best when sorting with the crowding metric. Evaluation of the remaining subjects in the front is performed by calculating the distance to its nearest neighbours. Advances in Production Engineering & Management 11(4) 2016 371 Klancnik, Hrelja, Balic, Brezocnik 1: Randomly generate and evaluate initial population of parents Po. 2: Prepare empty initial population of descendants Qo. 3: Setup t = 0. 4: Until stopping criterion has been fulfilled, repeat: 5: Setup t = t + 1. 6: Unite two old populations of parents and descendants: Rt-i = Pt-i % Qt i. 7: Perform non-dominated sorting on Rt i and determine fronts Fi, i = 1,2,... 8: Prepare new empty population of parents Pt. 9: Put into population Pt the first i front fronts still fitting whole into it. 10: With the use of crowding distance metric sort front F+i no more fitting whole into population Pt 11: Add the least crowded subjects from Fi+1 to population Pt 12: Generate population of descendants Qt from parent population Pt by using tournament selection, crossover and mutations. 13:_Evaluate subjects from population of descendants Qt._ Fig. 2 Pseudocode of the NSGA-II algorithm When optimising k criteria, k subjects are first sorted according to growing evaluating values f for each criterion j = 1, 2,., k. In the next step for each subject I, the distance between its two neighbours u and v is calculated according to the following equation: dj(i)=--(13) 7V y j?max _j?min where fjnax and f™111 are the maximum and the minimum values of the j-th criterion, respectively. The following must apply: fj(u) < fj(i) < fj(y) (14) As already mentioned, the highest possible distance is assigned to the two extreme subjects (with respect to the j criterion). For all remaining subjects, the crowding metric for subject i is equal to the sum of those distances according to all criteria: k c(i) = Yjdj(i) (15) 7 = 1 The parent population Pt is obtained in the manner described above. In the next steps, the descendant population is generated from that population by using tournament selection, crossover and mutation. Out of two random subjects, the subject classified into the front with lower consecutive value, wins in the selection. If two subjects taking part in selection come from an equal front, then the subject better evaluated by the crowding metric is chosen. The descendant population is designated Qt. Each subject from the descendant population is evaluated in the last step of the NSGA-II algorithm. That estimate is used in the next generation for non-dominated sorting of subjects. The algorithm is repeated until the maximum number of generations has been reached or until another stopping criterion has been reached. 3. Results and discussion 3.1 Modelling of turning process by the Gravitational Search Algorithm For processing experimental data the GSA method was used for the purposes of turning system modelling. For the modelling process the following control parameters of the GSA optimization algorithm were used: • Number of mass bodies: 210 • Number of iterations: 1000 • Size of mass body (dimension): 10 • Power of Euclidean distance between mass bodies: 1 372 Advances in Production Engineering & Management 11(4) 2016 Multi-objective optimization of the turning process using a Gravitational Search Algorithm (GSA) and NSGA-II approach • Power factor a: 20 • Gravitational constant Go: 100 In particular, the two random factors a and Go must be pointed out. Both of them are interrelated and form the refreshment of iteration within the scope of the gravitational constant, where the Go value is the initial value of the gravitational constant and is changed exponentially on the basis of factor a which, in combination with the ratio of consecutive iteration and total number of iterations, affects the contingency of gravitational constant values additionally during each iteration. The mass particle size and/or the properties of solutions written in the mass particle are represented with a dimension which also defines the number of coefficients in the polynomial representing the model of the individual machining operation. A prerequisite for successful execution of modelling was correct selection of the polynomial optimized by means of GSA. Eq. 16 shows the form of a ten-coefficient polynomial which proved to be a good choice for prediction of output variables. /(x1,x2,x3) = kt + k2 •x1 + /c3 • x2 ■ x3 •x-l • x2 + k6 ■ -x3 (16) + kj ■ x2 ^X3 + ■Xi ■Xi + ■ x2 ■ X2 + ^10 ■ ■ x3 where x1 is cutting speed Vc , x2 is feed rate f and x3 is depth of cut ap. The polynomial represents the core of the system round which the optimization process is designed. Input and output experimental values are inserted into the polynomial, while the optimization process optimises the coefficients in such a manner that the predictions given by the polynomial model are as near the experimental values as possible. Out of all the experimental measurements performed, fifteen measurements were chosen at random for the process of optimising of the polynomial coefficients, while the remaining 5 measurements were used for result verification. The results of optimization are polynomial models for the calculation of cutting force, surface roughness and tool wear (Eq. 17-19). Fc = 523.66 - 1.89 ■ x1 - 813.18 ■ x2 + 48.88 -x3 + 3.05- x1^x2 + 0.29 ■x1■x3 + 1439.13 ■ x2 ■x-j + 0.0014 ■ x1 ■x1 - 1257.73 ■ x2 ■X;, (17) - 43.70 ■x, l3 Ra = 3.74 - 0.01 • x1 -4.25-x2 + 1.53-x3 -0.01- x1 -x2 - 0.001 • x1 -x3 - 1.03 • x2 -x3 + 0.00002 • x1 •x1 + 68.81 -x2 -x2 -0.31-x3 -x3 T = 294.77 - 0.95 • x1 - 242.17 • x2 -4.60-x3 + 0.24- x1 -x2 + 0.02-xx -x3 - 9.20 • x2 -x3 + 0.0007 • x1 •x1 + 287.75 • x2 -x2 -5.37-x3 -x3 (18) (19) Testing of the developed models was recorded as a comparison of experimental values, predictions and percent deviation; experimental value ^ prediction ^ percent deviation: • Principal cutting force (Eq. 17) Lowest deviation with measured value: 550.848 N ^557.877 N ^1.28 % Highest deviation with measured value: 174.024 N ^189.579 N ^8.94 % Average deviation of all twenty measurements with ten predictions: 3.57 % • Surface roughness (Eq. 18) Lowest deviation with measured value: 2.32 [im ^2.31 [im ^0.43 % Highest deviation with measured value: 1.31 [im ^1.48 [im ^13.30 % Average deviation of all twenty measurements with ten predictions: 4.02 % • Tool life (Eq. 19) Lowest deviation with measured value: 10.93 min ^ 10.77 [im ^1.49 % Highest deviation with measured value: 28.43 min ^ 26.45 min ^ 6.97 % Average deviation of all twenty measurements with ten predictions: 4.50 % In this case, the deviation was within the acceptable ten percent. In any case, the lowest possible deviations are desirable, but the size of deviations depends primarily on the number of Advances in Production Engineering & Management 11(4) 2016 373 Klancnik, Hrelja, Balic, Brezocnik measurements. The experiment assured only 20 results of finishing machining, therefore, the final average of results can be assumed to be acceptable. 3.2 Multi-objective optimization by the NSGA-II algorithm On the basis of models created by means of the GSA algorithm and outlined in the previous subsection, in this paragraph, a set of optimal solutions for the chosen machining operation will be searched for. In our case, the resulting models for force, roughness and tool resistance represent three objective (criteria) functions having to satisfy also the limitations stated below. In this research, the aim of multi-objective optimization by the NSGA-II method was to determine such a combination of independent input variables (i.e. cutting force, feed rate and cut depth) at chosen intervals, that optimal values of criteria functions will be reached and the limitations prescribed satisfied. On the basis of theoretical calculations of limitations, and on the basis of experimental data, the optimization limitations of cutting force Fc, roughness Ra and tool resistance T were determined: • Fc < 450.0 N • 1.0 < Ra <1.6 |im • 15.0