Study of a novel cathode tool structure for improving heat removal in electrochemical micro-machining

Electrochimica Acta 75 (2012) 94–100

Contents lists available at SciVerse ScienceDirect

Electrochimica Acta journal homepage: www.elsevier.com/locate/electacta

Study of a novel cathode tool structure for improving heat removal in electrochemical micro-machining Jing Wu, Hong Wang ∗ , Xiang Chen, Ping Cheng, Guifu Ding, Xiaolin Zhao, Yi Huang Science and Technology on Nano/Micro Fabrication Laboratory, Research Institute of Micro/Nano Science and Technology, Shanghai Jiao Tong University, Shanghai 200240, People’s Republic of China

a r t i c l e

i n f o

Article history: Received 7 February 2012 Received in revised form 20 April 2012 Accepted 20 April 2012 Available online 4 May 2012 Keywords: EMM Hollow electrode Heat removal Multi-physics Pulse rate

a b s t r a c t In this paper, the effects of a hollow structure of cathode tool and radial flow of electrolyte on heat removal are studied using COMSOL software. A multi-physics (electrical field, flow field and heat transfer) model, which implement the bilateral interactions with customized PDE, is proposed to simulate the EMM process with a moving cathode tool. The ALE method is used to track the moving interface. The simulations show that the electrolyte flow rate and flushing time are the most important factors for the heat removal. A hollow structure of the cathode tool is proposed to improve the stability of flow rate, which will remove the heat generated during the EMM process effectively. © 2012 Elsevier Ltd. All rights reserved.

1. Introduction The growing demand for micro parts and molds has increased the importance of micro-machining technologies such as micro electrical discharge machining, micro electrochemical machining, and micro-electrochemical discharge machining. Among these processes, electrochemical micro machining (EMM) is used widely as it allows manufacturers to shape hard metals at a high material removal rate, without affecting the tensile strength of the workpiece material and its other physical properties, while ensuring a low surface roughness [1]. Although the process of EMM is difficult to predict, many mathematical approaches such as boundary elements method (BEM), finite difference method (FDM) and finite element method (FEM) have been used to analyze the electrical and fluid properties in the process. And many numerical models have also been proposed, ranging from the most commonly used potential models to multi-ion transport and reaction models [2,3]. Deconinck et al. [4] simulated the influence of the temperature on the uniformity or copying quality of the removal rate. Hourng et al. [5–7] calculated the temperature distributions in the electrolyte in one or two dimensions. Bhattacharyya and Munda [8] demonstrated that a significant metal removal rate and minimal overcutting could be obtained with the proper operating potential and electrolytic concentration and sufficient processing time. Kock

∗ Corresponding author. E-mail address: [email protected] (H. Wang). 0013-4686/$ – see front matter © 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.electacta.2012.04.078

et al. [9] verified the effects of pulse width on fabrication accuracy. Lee et al. [10] demonstrated the effect of pulse rate and the inter-electrode gap (IEG) value on the concentration distribution in IEG during the EMM process. It is well known that the EMM process is influenced by numerous factors such as temperature, IEG, electrolytic concentration, and operating potential. In the process, several physics are involved and flow field play an important role in the conduction of current between electrodes and maintaining the stability of the environment. The structure of cathode is a critical factor for the flow field excepting for the machining shape. The literatures mentioned above all used a solid cathode tool in the process. However, these solid tools have several disadvantages such as low copying accuracy and asymmetry. A hollow structure of the cathode tool, resulting in high copying curacy and enhanced symmetry, has been adopted in traditional large scale machining. Since the adverse influence of decrease in flow rate is more serious in EMM, the adoption of hollow cathode tool would be a significant improvement. In the following sections, the paper discusses a multi-physics model that presents some comparative simulations to study their sensitive relevance on temperature distribution in EMM process. 2. Theoretical model EMM is a current controlled electrolytic process. Typically, a workpiece (anode) is connected to a tool (cathode) and immersed in an ionic solution in a reactor. When a potential difference is applied between electrodes, the oxidation and reduction reaction

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

95

will take place on anode and cathode, respectively. The metallic ions of the anodic dissolution will migrate toward the cathode and form the metallic hydride and oxide after reacting with the OH− . At the same time, oxygen is generated on the anode while hydrogen is generated on the cathode. By feeding the shaped cathode toward the workpiece, the reverse shape of the tool is copied. In order to attain a more specialized study of the temperature variation in the EMM process, the aspects of concentration and gas evolution are not considered here. As the electrode reaction and current density are strongly dependent on temperature, and heat, generated by the electrochemical reaction, is mainly removed by the electrolyte flow. It was reasonable to consider the electrolyte only in the view of obtaining the temperature distribution. Electrodes are usually good heat conductors, approximately 100 times better than the electrolyte and it is logical to regard the electrode–electrolyte interface as outer boundaries (reference temperature). During this process, several controllable parameters such as, operating potential, flow rate of electrolyte, IEG, and electrolyte composition, play an important role on the copying quality. The structure of cathode tool, which is a critical factor for maintaining flow rate in the process, is transformed to hollow from solid. 2.1. Model assumptions The EMM process involves several interactive physics, such as flow field, electric field and heat transfer. First, a model of electrodes with given dimensions is needed to attain the potential distribution. And then an electrolyte model with potential U, temperature T and electrolyte velocity v is developed according to the experiments. An established heat transfer model acts as a bridge coupling the two mentioned above. In order to obtain coupling physics results, a multi-domain method is used. To simplify the proposed model, several assumptions are made as follows: (1) The EMM process is ideal. The anodic dissolution will take place as long as a potential difference is applied between electrodes, and it depends only on temperature. (2) The model is isothermal. The thermal conductivity of the electrolyte was assumed to be a constant.

Fig. 1. Geometry model with a downward feeding cathode tool; the electrolyte flow region is noted as gray area: (a) initial model and (b) simplified simulation model.

The anode interface and cathode interface are taken as equipotential surfaces, denoted as Va and Vc , for solving the charge conservation

∇ · (∇ U) = 0

Va , Vc are imposed as constant potential (V). Meanwhile, a temperature dependent model, in which the electrical conductivity is linearly related to temperature T, is selected to obtain the potential distribution. The relation between  and T is introduced:  = [1 + ˛(T − Tref )]

2.2. Geometry model The geometry used in this study is presented in Fig. 1, a 2D channel, comprised a cathode surface (no. 3, 4, 5, 7, 8 and 9), an anode boundary [12], containing the electrolyte flowing between them. The channel is 150 ␮m in width and has a length of 600 ␮m. The inter-electrode gap (IEG), is set to 20 ␮m. Radius of the hollow structure is 50 ␮m. In addition, the hollow tool feeds down with a constant rate vc in the process. The anode surface, which is in contact with the electrolyte, is fully active with a constant potential applied to it; while the cathode tool is side-insulated, i.e., boundaries 3 and 9 are electrically insulated from the electrolyte. The electrolyte flows in radial flow at a uniform rate profile, i.e., flows in from AB and exits from CD and EF. In order to avoid the convergence in the solving process, the Heaviside function [10] of the velocity distribution is set as the initial value at the inlet, AB.

(2)

(3)

where Tref is the reference temperature and a conductance constant ˛ indicates the degree of temperature dependence. A typical potential distribution is obtained, shown in Fig. 2. In order to decrease the influence of stray current density in processing, the sidewall insulation is adopted.

2.3. Electrical model The potential between the electrodes comply with Laplace equation,

∇2 = 0 i.e., the potential field is a steady current field.

(1)

Fig. 2. Potential (V) distribution in the vicinity of the cathode tool after 1 s of metal dissolution processing; colored curves represent isopotential curves. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

96

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

Fig. 3. Electrolyte flow field in different geometry structure in the processing. The inlet velocity of is 0.5 m/s .The color-scales represent the magnitude of the velocity: (a) the solid structure cathode tool, side flow and (b) the hollow structure cathode tool, radial flow. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

In electrical model discussed above, the heat released in the EMM process, which will be discussed in Section 2.5 in detail, is bound to heat the electrolyte up rapidly. In order to promptly remove the heat and prevent the electrolyte from boiling, an appropriate and effective flow field becomes a necessity.

left side term in equation accounts for heat removed by the convection effect. PJ (W/m3 ) represents Joule heating in the bulk of the electrolyte and it can be attained by Joule’s law, P J = E · j =  (∇ U)2 

2.4. Fluid flow model The electrolyte velocity v is obtained by solving the incompressible laminar Navier-Stokes equations,

∇ ·v=0





∂v + v · ∇v ∂t

(4)

 = −∇ p + v

(5)

where p is the pressure (Pa),  is the electrolyte density (kg/m3 ) and  is the dynamic viscosity (kg/(m s)). The effects of different structure of cathode tool (solid and hollow) on flow rate are simulated, shown in Fig. 3. The flow field using solid cathode tool is presented in Fig. 3(a). The workpiece (anode) is connected to the solid tool (cathode) and immersed in an ionic solution in an electrolyzer. The electrolyte is pumped from the left side, in other words, the type of flow is side flow. The electrolyte flow passage between electrodes gradually becomes a small ‘U’ channel from a, relatively, wide rectangle channel as the tool feeding. Velocity declines rapidly to 0.1 m/s, accounting for one fifth of the initial value, which would seriously impact the convection effect on heat removal. Hence, it is necessary to modify the solid structure to maintain the stability of the flow rate. As shown in Fig. 3(b), a hollow structure, allowing radial flow, is adopted in this simulation to maintain the stability of flow rate. It is known that dynamic viscosity, density and many other flow properties vary with temperature. And a mathematical relation needs to be developed for studying the interaction between these properties and temperature T. In our simulations, values of these properties at different temperatures are fitted using MATLAB. The equations are obtained as follows:  = 0.9498 · e−0.0524·T + 0.08352e−0.01115·T  = 1101 · e

−0.001076·T

+ (−101)e

−0.01208·T

(6) (7)

2.5. Electrolyte heat transfer model In electrolyte, temperature T is found by solving internal energy balance, Cp v · ∇ T = ∇ · (k∇ T ) + P J + P dl

(8)

where Cp is the specific heat capacity (J/(kg K)), k is the thermal conductivity and v is the velocity of the electrolyte flow (m/s). The

(9) conductivity(S m−1 ),

where is the electrical U is the potential difference applied over the electrodes and E is the electric field. In addition, there is another heat source Pdl (W/m3 ), which represents the heat generated by the electrochemical reaction in the double layers (DLs). Maximum heat is generated when enthalpy changes are neglected, yielding P dl = j

(10)

where  is the overpotential (V) and j is the reaction current density (A/m2 ). This polarization behavior is liberalized with respect to overpotential for electrode reactions: j=

(V − U) − Q R

(11)

where R is the polarization resistance ( /m2 ); Q is the onset of polarization (V), the minimum overpotential to initiate the reaction. From the correlation between the reaction current density j and the temperature T [22], it is reasonable to introduce the following equation: j = j[1 + ˇ(T − Tref )]

(12)

where ˇ is a constant that determines the magnitude of temperature dependence, j is the instant reaction current density varying with T. The modification of reaction current density j ultimately reflects on the heat production from the reaction, and the heat boundary condition transfer to the expression as follows: H  = H[1 + ˇ(T − Tref )]

(13)

where H , and H represent the instant heat production and the initial value, respectively. From measurements [4], it is known that H approximately equals 106 (W/m2 ), which is set as the initial value of H to simplify the numerical model. In the EMM process, the heat removal is controlled by the convection flush of electrolyte flow and by heat conduction resulting from temperature difference. However, the contribution of heat conduction is less than 0.1% [4], and the interfaces between the electrodes and electrolyte are considered as outer boundaries (temperature reference). In this sense, it is logical to consider the convection influence only in the view of obtaining the temperature distribution during the EMM process. A typical simulation result for the coupling of the temperature and flow field can be obtained, shown in Fig. 4. As discussed in Sections 2.3 and 2.4, the velocity at the inlet is uniformly set to a

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

97

Fig. 4. Temperature (◦ C) distribution in the vicinity of the cathode tool after 1 s of metal dissolution processing.

constant vertical vin , and no-slip conditions are imposed on other boundaries. A constant pressure is applied on the outlets. The radius of the hollow structure and the inlet velocity are set as 50 ␮m and 0.5 m/s; other parameters used can be found in Section 2.7.

In order to maintain a relatively stable IEG, the cathode tool moves at a uniform velocity, which is chosen to be 1.2 × 10−5 (m/s) according to the experiments (see Table 2). 3. Simulation results and discussion

2.6. Moving mesh (ALE) In order to track the movement of the anode boundary and the feed of the cathode tool, the computational geometry has to be realigned with the new position of the electrodes timely. Meanwhile, the computational mesh alters inevitably. And the arbitrary Lagrangian–Eulerian (ALE) method, which is an intermediate between the Lagrangian and Eulerian methods, is selected for implementing the deformed mesh and moving boundaries, as it allows moving boundaries without the need for the movement of mesh. The anode boundary moves downward at the rate of va , which can be obtained by the following equation:

va = ωj

(14)

Several comparative simulations were performed with the multi-physical model shown in Fig. 1. The simulations were run in transient mode using the Heaviside function of velocity to avoid the convergence in solving. 3.1. Influence of hollow structure The influence of the hollow structure is investigated. The comparative results are shown in Fig. 5. When a solid tool or a small radius of a hollow tool was chosen, the weak convection cannot effectively remove the heat in IEG during the process. The black and red curves in Fig. 6 show that in the above two cases, the temperature rises sharply and exceeds 35 ◦ C within 600 ␮s, exhibiting an ascending tendency throughout. T increases slightly and presents

(A/m2 ),

where j is the current density ω is the volume equivalent (m3 /(A s)), and  is a dimensionless constant that represents dissolution efficiency. The feeding rate of the cathode tool tip vc , and the workpiece dissolution rate va are imposed as Dirichlet boundary conditions to implement ALE method. 2.7. Parameters and boundary conditions The parameters used in simulations are listed in Table 1. The dynamic viscosity and density are considered as temperature dependent variables, as mentioned in Section 2.4. The inlet velocity is 0.5 m/s, and the type of flow is radial. The outer reference temperature and the initial temperature are assumed to be 20 ◦ C (room temperature). The temperature dependence factor ˛ for electrical conductivity of the electrolyte is taken as 0.02, which is usually adopted for electrolyte solutions in EMM [20]. Table 1 Parameters and coefficients of the numerical model. Parameters

Values

Potential difference (˚/V) Pulse period (␮s) IEG value (␮m) Electrolyte electric conductivity (/S m−1 ) Thermal conductivity (/W m−1 K−1 ) Radius of hollow cathode (radi/␮m) Duty cycle ( /[1]) Temperature reference (Tref /◦ C) Inlet flow rate of electrolyte (vin /m/s)

5 200 20 1.542 (initial) 0.58 10, 30, 50, 70 0.2, 0.5, 0.8 20 0.001, 0.01, 0.1, 0.5

Fig. 5. Temperature distribution with different radius of hollow cathode tools after 1 s of processing; uin = 0.5 m/s: (a) solid cathode tool, side flow; (b) radius = 10 ␮m, radial flow; (c) radius = 30 ␮m, radial flow; and (d) radius = 70 ␮m, radial flow.

98

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

Table 2 Boundary conditions of the model. Boundary settings of the electrical model (ec) Boundary no. 1, 2, 3, 9, 10, 11, Boundary cond. Electric insulation Setting n·J=0

4, 5, 7, 8 Electric potential V = V0

12 Ground V=0

Boundary settings of the fluid model (spf) Boundary no. Boundary cond. Setting

1, 11 Outlet p=0

2, 3, 4, 5, 7, 8, 9, 10 No slip wall v=0

1, 2, 3, 9, 10, 11, 12 Thermal insulation −n · (− k · ∇ T) = 0

6 Temperature reference Tref

6 Inlet v = −vin · n

Boundary settings of the heat transfer model (ht) Boundary no. 4, 5, 7, 8 Boundary cond. Heat flux Setting H = H[1 + ˇ(T − Tref )]

flat trends only when the radius of hollow structure increases. One can also see that T goes up mildly in the beginning and then becomes stable at 25 ◦ C when the radius increases to 30 ␮m or higher, which can provide a stable temperature environment for the process. Although the convection effect can be enhanced gradually with an increase in the radius of the hollow structure, the adverse results are much more conspicuous after the radius exceeding 30 ␮m, as shown in Figs. 6 and 7. An increase in radius over 30 ␮m brings a decrease of one or tenths in temperature only; while the bulge in anode profile becomes much larger, even reach 10 ␮m at a radius of 70 ␮m. From above, one can conclude that the hollow structure improves the stability of flow rate and the radius cannot be too large to prevent the seriously adverse effects of bulge. 3.2. Influence of inlet velocity Fig. 8 investigates the influence of the inlet velocity vin on temperature variation during the process. From the results, one can see that T rises rapidly, as large as 10 ◦ C only in 200 ␮s, which is due to the weak convection. And the growth rate of temperature even shows a rising trend over the processing time. In contrast, the growth decreases a lot when the inlet velocity vin was imposed as 0.5 m/s. And the rising trend also becomes flat, which can facilitate a stable environment for improving the copying quality. The dotted curves represent the decrease in T from 25 ◦ C, and there is almost no decline in T if the flow rate is under 0.1 m/s. When the velocity of 0.5 m/s is applied, the time for the

Fig. 6. Temperature variation curves with different radius of hollow cathode tools at a certain site in IEG over processing time; uin = 0.5 m/s.

Fig. 7. Anode profile opposite the hollow parts of cathode tools after 1 s of processing; uin = 0.5 m/s.

temperature decreasing from 25 ◦ C to 21 ◦ C is about 200 ␮s, which can be seen as a reference for selecting the flushing time. And then

Fig. 8. Temperature rise according to different flow rate over processing time at a certain site in IEG: the different colors correspond to different flow rates; solid curves represent temperature increase; and dotted curves represent temperature decrease. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

99

4. Conclusions

Fig. 9. Variation in temperature at a certain site in IEG over processing time under different duty cycles; pulse rate = 200 ␮s: solid curves represent uin = 0.5 m/s; dotted curves represent uin = 0.1 m/s.

the introduction of pulse rate becomes necessary to prevent the excessive increase in T in the process.

3.3. Influence of duty cycle

Acknowledgments

Duty cycle (duty factor) is the ratio of on/off time, as shown in equation: Duty cycle =

ton ton + toff

In this work, a novel structure of cathode tool using in EMM was presented and a simplified and specialized model was proposed. The laws involved in the coupled model, such as charge conversation, internal energy balance and the Navier-Stokes equations for incompressible viscous flow, are solved using a multi-domain approach, which implement the interactions by partial differential equations (PDE). The ALE method is applied to describe the moving electrodes interface and the quality of the grid in transient is also ensured by the automatically re-mesh technology. The simulations confirm that electrolyte flow rate and convection time are of paramount importance in heat removal during the EMM process, as most of the generated heat is removed by electrolyte flow. It is also shown that the hollow cathode with an interior channel, compared with the solid cathode tool, is more conductive to fluid flow. The difference of the damping effect of the fluid flow in the gradually shrinking passage with the two cathode structures clearly illustrates the necessity to replace the solid cathode with a hollow one. Although most of the heat in IEG can be removed, the temperature will rise due to the accumulation of residual heat. And the convection would require abundant flushing time to confine the residual heat to a certain range. Pulse rate plays an important role in providing sufficient heat removal time. Pulse off time, which must be longer than the minimal time for bringing down the temperature, is found an critical factor to maintain plenty of flushing time.

(15)

The importance of pulse rate applied to operating potential is investigated. Variations in the temperature dependent parameters are acceptable and can be regarded as a constant when temperature T varies in a range of 5 ◦ C. In this study, the processing time for an increase of 5 ◦ C in temperature is adopted as the maximal pulse width, and the time for T dropping to Tref from 25 ◦ C, is selected as the minimum pulse off time. The two eventually determine the value of the duty cycle. The results (mentioned in Section 3.2) suggest that the temperature T rise 5 ◦ C after 100 ␮s of EMM processing (shown by solid curves) and can be cooled down to 21 ◦ C within 100 ␮s (shown by dotted curves).Therefore, it is reasonable to set the pulse period 200 ␮s in the simulation. The optimization of the duty cycle is also investigated. As shown in Fig. 9, in contrast to an increase of 17 ◦ C at the flow rate of 0.1 m/s, temperature T only rises approximately 6 ◦ C when vin = 0.5 m/s (duty cycle = 0.8), which could mainly maintain the stability of the electrolyte properties. This result verifies the necessity of the adoption of hollow structure to avoid significant decrease of flow rate in the EMM process. In addition, T has a big rise and is trapped in IEG when duty cycle exceeds 0.5. The difference in T with different duty cycle (0.8, 0.5) reached 11 ◦ C within only 1000 ␮s of the processing time. And it probably increases to a value that cannot be ignored as the processing time extends. So it is important to confine the duty cycle to prevent the excessive increase in T. However, undesirable decrease of dissolution rate would occur if the pulse-on time were too short. In order to balance the requirements of the two, medium pulse rate and duty cycle become fairly necessary.

This research was supported by the Advanced Research Project under Grant (625010105), the Program of Science and Technology Commission of Shanghai Municipality under Grant (11DZ2290203) and the Cultivating Foundation for Innovation Engineering Program of the Ministry of Education of China under Grant (708037). References [1] D.K. Chung, H.S. Shin1, M.S. Park, B.H. Kim, C.N. Chu, International Journal of Precision Engineering and Manufacturing 12 (2011) 371. [2] D. Deconinck, S.V. Damme, J. Deconinck, Electrochimica Acta 60 (2012) 321. [3] D. Deconinck, S.V. Damme, J. Deconinck, Electrochimica Acta 69 (2012) 120. [4] D. Deconinck, S.V. Damme, C. Albu, L. Hotoiu, J. Deconinck, Electrochimica Acta 56 (2011) 5642. [5] L.W. Hourng, C.S. Chang, Journal of Applied Electrochemistry 24 (1994) 1170. [6] M.J. Noot, R.M.M. Mattheij, Mathematical and Computer Modelling 31 (1997) 77. [7] C. Sun, D. Zhu, Z. Li, L. Wang, Finite Elements in Analysis and Design 43 (2006) 168. [8] B. Bhattacharyya, J. Munda, International Journal of Machine Tools & Manufacture 43 (2003) 1301. [9] M. Kock, V. Kirchner, R. Schuster, Electrochimica Acta 48 (2003) 3213. [10] Y.M. Lee, S.J. Lee, C.Y. Lee, D.Y. Chang, Journal of Power Sources 193 (2009) 227. [12] J.A. McGeough, Principles of Electrochemical Machining, Chapman and Hall, London, 1974. [20] J. Newman, K.E. Thomas-Alyea, Electrochemical Systems, 3rd ed., John Wiley & Sons, Inc., Hoboken, NJ, 2004. [22] K.A. Connors, Chemical Kinetics. The Study of Reaction Rates in Solution, VCH Publishers, 1990.

Further reading [11] M. Purcar, Development, Evaluation of Numerical Models and Methods for Electrochemical Machining and Electroforming Applications, Universiteit Brussel, 2005. [13] J. Jacquemin, P. Husson, A.A.H. Padua, V. Majer, Greem Chemistry 8 (2006) 172. [14] T. Wuilbaut, Algorithmic Developments for a Multiphysics Framework, Ph.D. Thesis, Université Libre de Bruxelles, von Karman Institute for Fluid Dynamics, 2008. [15] D. Zhu, H.Y. Xu, Journal of Materials Processing Technology 129 (2002) 15. [16] H.S. Shin, B.H.K. Chong, N. Chu, Journal of Micromechanics and Microengineering 18 (2008) 075009.

100

J. Wu et al. / Electrochimica Acta 75 (2012) 94–100

[17] R. van Tijum, P.T. Pajak, Simulation of production process using the multiphysics approach: the electrochemical machining process, in: Excerpt from the Proceedings of the COMSOL Conference, Hannover, 2008. [18] J. Kozak, Journal of Materials Processing Technology 76 (1998) 170. [19] M.A.H. Mithu, G. Fantoni, J. Ciampi, International Journal of Advanced Manufacturing Technology 55 (2011) 921.

[21] G.H. Liu, Y. Li, X.P. Chen, S.J. Lv, Advanced Materials Research 60 (2009) 380. [23] E.S. Lee, T.H. Shin, International Journal of Precision Engineering and Manufacturing 11 (2010) 113. [24] R. Schuster, V. Kirchner, P. Allongue, G. Ertl, Science 289 (2000) 98.