- Email: [email protected]
Proposal of a new design for valveless micropumps
Scientia Iranica B (2011) 18 (6), 1261–1266
Sharif University of Technology Scientia Iranica Transactions B: Mechanical Engineering www.sciencedirect.com
Proposal of a new design for valveless micropumps H. Afrasiab, M.R. Movahhedy ∗ , A. Assempour Center of Excellence in Design, Robotics and Automation, Department of Mechanical Engineering, Sharif University of Technology, Tehran, P.O. Box 11155-9567, Iran Received 7 February 2011; revised 9 August 2011; accepted 8 October 2011
KEYWORDS Valveless micropump; Finite element simulation; Fluid–Structure Interaction (FSI) analysis; Arbitrary Lagrangian–Eulerian (ALE) method.
Abstract A new design for a valveless micropumping device has been proposed that integrates two existing pumping technologies, namely, the wall induced traveling wave and the obstacle-type valveless micropump. The liquid in the microchannel is transported by generating a traveling wave on the channel, while the placing of two asymmetric trapezoid obstacles, along the centerline of the channel inlet and outlet, leads to a significant (up to seven times) increase of the net flow rate of the device. The effectiveness of this innovative design has been proved through a verified three-dimensional finite element model. Fluid–Structure Interaction (FSI) analysis is performed in the framework of an Arbitrary Lagrangian–Eulerian (ALE) method. © 2012 Sharif University of Technology. Production and hosting by Elsevier B.V. All rights reserved.
1. Introduction The growing requirement of controlled fluid transfer in many MEMS and BioMEMS applications has made the development of microfluidic systems an area of increasing interest in recent years. Dispensing drugs and other therapeutic agents, fuel delivery in miniaturized fuel cells and provision of refrigeration liquid in the cooling of microelectronic systems are some examples of these applications. Micropumps are essential components of microfluidic systems. Several different types of micropump have been developed up to now, valveless piezoelectric micropumps of which are in wide practical use amongst which valveless piezoelectric micropumps are in wide practical use, due to their simple structures, high flow rate output and ability to conduct particles such as red blood cells, polymers and proteins. These kinds of micropump are driven by a piezoelectric element bounded to a flexible membrane, and are equipped with a
∗
Corresponding author. E-mail address: [email protected] (M.R. Movahhedy).
1026-3098 © 2012 Sharif University of Technology. Production and hosting by Elsevier B.V. All rights reserved. Peer review under responsibility of Sharif University of Technology. doi:10.1016/j.scient.2011.11.023
pressure chamber for generating a pressure difference between the inlet and outlet, and no interior moving mechanical part is present in them. Instead, flow directing components like nozzle/diffuser elements or asymmetric obstacles are utilized to control rectification and obtain the net volume of the fluid flow [1,2]. Though piezoelectric micropumps can transport liquid at a high flow rate and with excellent controllability, the presence of the pressure chamber makes further miniaturization difficult and sometimes impractical [3]. Seeking a remedy for this problem, recently a mechanical micropumping device was fabricated and presented in [3], which uses a vibrating microchannel wall to drive the fluid flow while no pressure chamber is present. As reported in [3], on-generating a traveling wave on the channel wall, clear liquid flow could be induced within the microchannel. The transportation of the liquid by a traveling wave is controlled by changing the magnitude and frequency of the generated wave. A range of frequencies between 50 and 600 kHz were considered in [3] for oscillation of the microchannel wall for a maximum oscillation magnitude of 142.6 nm. A net outlet flow rate of 1.20 nL/s for the vibration frequency of 100 kHz was reported in [3] for a microchannel of 200 µm width and 300 µm height. The simple structure of this Traveling Wave Micropump (TWM) makes it advantageous for many microfluidic systems. In this paper, we propose a simple geometric modification to the traveling wave micropump so as to increase its net outlet flow rate. This modification is inspired by the action mechanism of the obstacle-type valveless micropump. The paper is organized as follows. The working principle of the new micropumping device is described in Section 2. The finite element modeling procedure is presented in Section 3. In Section 4, the results of FEM simulations are presented
1262
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
Figure 1: Working principle of the obstacle-type valveless micropump [5].
micropump. Now, we place two trapezoid obstacles along the centerline of the microchannel, one in the channel inlet downstream of the electrodes array and the other in the channel outlet upstream of the electrodes array, as shown in Figure 2. These obstacles form two equivalent planar nozzle/diffuser rectifying elements, similar to those described for the obstacletype valveless micropumps above. Consequently, a directiondependent flow resistance is provided, and an increase in the net flow rate of the microchannel is expected based on the previous experience with obstacle-type micropumps. In the next sections, we examine the effect of the proposed geometric modification on the net outlet flow rate of the microchannel, using a finite element simulation of both with- and withoutobstacle microchannels. 3. Finite element method simulation
Figure 2: A schematic of the proposed micropump.
and discussed. Finally, the concluding remarks are reviewed in Section 5. 2. Working principle of the new micropump In order to shed light on the action mechanism of the proposed micropump, first we briefly review the working principle of the obstacle-type valveless micropump presented e.g. in [4,5]. This kind of micropump is composed of a pump chamber with an oscillating diaphragm and two nozzle/diffuser flow rectifying elements, one directed from the inlet to the pump chamber, and the other from the pump chamber to the outlet, as illustrated in Figure 1. Each element is formed by two equivalent regions between the symmetrical trapezoid obstacle and the side-walls of the microchannel. Actuating the diaphragm causes a continual periodic increase and decrease in the volume of the pumping chamber (Figure 1). When the volume of the pumping chamber increases, the pressure in the chamber decreases, and more fluid enters through the nozzle-diffuser element on the left (inlet) relative to that on the right (outlet). This is because the element on the left acts as a diffuser, which poses less flow resistance than the nozzle on the right. Conversely, when the volume of the pumping chamber decreases, more fluid exits through the element on the right, which now acts as a diffuser. This results in a net pumping action from left to right, in Figure 1, as the diaphragm vibrates up and down. Now, we turn our attention to the new micropump. A schematic of this device is illustrated in Figure 2. In this device, the fluid flow is driven by the vibrating microchannel wall and actuated by a piezoelectric layer. The top ceiling wall of the microchannel is composed of a silicon wafer whose surface is covered by a piezoelectric PZT thin film. Nine separate top electrodes are attached to the PZT film. The traveling wave is induced on the ceiling wall by applying sinusoidal voltages, with a phase difference of 2π /3, to the adjacent top electrodes. When a channel wall oscillates in the form of traveling waves, peristaltic motion is induced in the liquid beneath the microchannel wall [3]. After a period of oscillation, the fluid moves slightly forward from the initial position, due to its viscosity, and the net liquid can be transported by repeating the period of motion. Up to this point, the pumping principle of the device is exactly similar to that of the traveling wave
In order to investigate the effect of obstacle addition on the net flow rate of the micropump, numerical simulations of both with- and without-obstacle microchannels were carried out in this section, using the finite element method. As stated earlier, the fluid flow is driven by the vibrating microchannel wall. Since the fluid also plays a role in resisting the wall vibration, the solid wall vibration and fluid flow are coupled, and a full Fluid–Structure Interaction (FSI) analysis is required in order to investigate the performance characteristics of these devices [6,7]. For this purpose, the FSI analysis is performed, using an Arbitrary Lagrangian–Eulerian (ALE) approach in the framework of the finite element method. The simulations are performed by means of a verified in-house developed FSIALE code, which has been used previously, e.g. in [8,9]. The finite element modeling details are described in the following subsections. 3.1. The fluid flow governing equations The law of conservation of momentum and the continuity equation for a viscous Newtonian incompressible flow are Navier–Stokes equations. In the ALE approach to fluid–structure interaction problems, these equations are written in the ALE description to account for deformation of the fluid mesh in the FSI interface. These equations, along with appropriate initial and boundary conditions, are seen in the ALE description (see e.g. [10])
ρf
∂ v − 2µf ∇.∇ s v + ρ v − vm .∇ v + ∇ p ∂t χ = ρ bf
in Ω f × (0, T ) ,
(1a)
∇.v = 0 in Ω f × (0, T ) ,
(1b)
v (x, 0) = v0
(1c)
on Ω × {0} , f
v (x, t ) = vD (x, t )
on
f ΓD
× (0, T ) ,
n .σ = −pn + 2µf n .∇ v = t f
f
f
f
s
f
on
(1d) f ΓN
× (0, T ) ,
(1e)
where ρf , µf , v and p represent the fluid density, viscosity, velocity and pressure, respectively. vm is the mesh velocity and bf denotes the body force vector. The fluid domain is specified by Ω f , while the time interval of interest is between (0, T ). The ∂(·)
time derivative in the ALE frame is represented by ∂ t . The χ
symmetric tensor, ∇ s v = ∇ v + (∇ v)T /2, is called the rate of deformation (or strain rate) tensor. nf is the exterior normal to the fluid boundary, Γ f , and v0 is the initial condition for
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
1263
f
the velocity field. vD denotes the prescribed velocity on the ΓD portion of the boundary, and tf is the boundary traction on the f complementary portion, ΓN . Since the working frequency of the solid wall is very high (between 50 and 600 kHz), the time step size required for the finite element analysis of these micropumps is very small. Thus, in order to avoid the so-called small time step instability in the finite element analysis of the fluid part, we used a residualbased variational multiscale formulation for incompressible Navier–Stokes equations, in combination with stabilization parameters, as recommended in [11], for weak formulation of the fluid flow. 3.2. The solid governing equations In a standard Lagrangian description, the law of conservation of linear momentum for a solid continuum may be expressed with respect to spatial coordinate, x, as:
∂ 2u = ∇.σs + ρs bs in Ω s × (0, T ) , ∂t2 u (x, 0) = u0 (x) in Ω s × {0} , ∂u (x, 0) = u˙ 0 (x) in Ω s × {0} , ∂t u (x, t ) = uD (x, t ) in ΓDs × (0, T ) , ρs
(2a)
Figure 3: The sequential fluid–structure coupling algorithm.
ρf vm .∇ v takes the convection of mesh moving points into
(2d)
account. The solid structure is formulated and solved in the Lagrangian description, as indicated earlier in this paper. Different techniques of ALE mesh moving are presented in the literature. Here, we use an approach that is based on solving the linear elasticity equations [12]. In this approach, the equation governing the displacement of fluid mesh nodes can be written as:
σ .n = t on × (0, T ) , (2e) where ρs is the solid density, and u represents its displacement
∇.σ = 0 in Ω f × (0, T ) , (5) where σ is the Cauchy stress tensor. For each boundary, i, a
s
s
s
(2b) (2c)
ΓNs
field, whereas the body forces are given by vector bs . The symmetric second order tensor, σs , denotes the Cauchy stress tensor, and the solid domain is represented by Ω s . uD denotes the prescribed velocity on the ΓDs portion of the boundary, and ts is the boundary traction on the complementary portion, ΓNs . The outer unit normal to the solid boundary is specified by vector ns . Since the solid membrane deformation is small in most microfluid manipulating devices [8], we assume that the structural part of the FSI problem is governed by linear elasticity. For a linear elastic solid, we have:
σs = λs (∇.u) I + 2µs ε (u) , (3) T in which ε (u) = ∇ u + (∇ u) /2 is the strain tensor, and µs and λs are the Lamé coefficients. 3.3. Fluid–structure interface conditions Fluid–solid interface conditions consist of kinematic and dynamic constraints specified as follows on the FSI interface, ∂ Ω fs :
∂ us = vf ∂t
on ∂ Ω fs × (0, T ) ,
σs .n + σf .n = 0 on ∂ Ω fs × (0, T ) ,
(4a) (4b)
with n being the outer normal at the solid boundary. 3.4. The mesh motion algorithm In the Arbitrary Lagrangian–Eulerian (ALE) approach to fluid–structure interaction problems, two key ingredients are needed. The first is a technique to move the fluid mesh, so that it can track structure motion at the FSI interface. This deformation is extended into the fluid interior field, in order to avoid excessive distortion of the fluid mesh at the interface. The second ingredient is a flow formulation written in the ALE framework, as presented in Eq. (1). The time derivative in Eq. (1a) is taken within the ALE description and the term,
Dirichlet boundary condition: ui = ubi ,
(6) b
may be given. u is the boundary displacement vector that is either given a priori or computed by the solid structure equations. 3.5. Fluid–structure coupling method In this work, a partitioned strong coupling approach is employed that uses a sequential iterative algorithm, as illustrated in Figure 3. This iterative scheme for the coupled system of fluid and structure is repeated at each time step until convergence is reached, i.e. after the norms of relative change in the field variables between two consecutive iterations are smaller than the specified convergence tolerance. 3.6. Boundary and initial conditions The boundary conditions are as follows: Along the lateral and bottom walls of the microchannel and obstacles (fixed boundaries), the usual ‘‘no slip’’ boundary condition is prescribed for the fluid. Along the top wall of the microchannel, which is also the bottom of the microchannel ceiling (FSI boundary), the fluid–structure interface conditions given in Sections 3 and 4 are applied. The motion of the microchannel ceiling is prescribed on its upper face. The initial condition is zero flow velocity and structure (microchannel ceiling) displacement. 3.7. Verification of the numerical approach In order to verify the accuracy of the finite element approach, a portion of the without-obstacle microchannel (containing all electrodes, as well as 1 mm downstream and 1 mm upstream of the electrodes array) was simulated in three dimensions, for 20 complete periods of oscillation and for different values of oscillation frequency. 14 040 and 3120 tri-quadratic elements
1264
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266 Table 1: Without-obstacle microchannel geometry. Silicon wafer thickness (µm) Piezoelectric film thickness (µm) Distance between single electrodes centers (µm) Excitation voltage (V) Maximum oscillation amplitude of electrodes for frequency of 100 kHz (nm) Microchannel width (µm) Microchannel height (µm) Microchannel length (mm)
Table 3: With-obstacle microchannel geometry. 6.0 2.5 400 20 142.6
Distance between two obstacles (L) (mm) Obstacles’ length (L1 ) (µm) Microchannel width (W ) (µm) Obstacles’ width (W1 ) (µm) Obstacles’ divergence angle (α )
3.2 500 200 25 5°
200 300 5.2
Table 2: Solid physical properties. Material Silicon Piezoelectric
Density (kg/m3 ) 2300 7800
Elastic modulus (GPa) 170 64
Poisson’s ratio
Figure 5: A schematic top view of the micropump geometry.
0.215 0.31
Figure 4: Frequency dependence of the mean flow velocity.
were used to discretize the fluid and solid parts, respectively. We used 80 time steps per period of oscillation in all these simulations. The without-obstacle microchannel geometry is presented in Table 1. The working fluid is water with density of 1000 kg/m3 and viscosity of 0.001 kg/(ms). The solid physical properties are listed in Table 2. In each case, the mean flow velocity was calculated in the microchannel outlet. Figure 4 compares the results of these simulations with the experimental data reported in [3]. It can be seen that the finite element results are in good agreement with experimental data, which demonstrates the capability of the employed numerical model. 4. Results and discussion Once the numerical model is verified, it can be used for studying the performance characteristics of the proposed micropump. To this end, in this section, a simulation similar to that of Section 3 was carried out for the with-obstacle micropump, for wall oscillation frequency of 100 kHz. 16 448 and 4112 tri-quadratic elements were used to discretize the fluid and solid parts, respectively. A schematic top view of the micropump geometry is indicated in Figure 5, and detailed dimensions of the obstacles and microchannel are listed in Table 3. The dimensions of the obstacles are adopted from those presented in [4]. Other geometric dimensions of the withobstacle micropump that are not given in this table are the same as those of the without-obstacle micropump presented in Table 1. It is to be noted that due to the presence of the
Figure 6: Outlet flow rate variation for without- and with-obstacle micropumps.
obstacles, more energy is needed for suction and discharge of fluid into and from the micropump. Consequently, in the new micropump, a higher voltage is needed for deforming the microchannel wall that controls the fluid flow. An electrofluid–structural interaction simulation (see e.g. [13,14]) is needed for determination of the required increase in the operating voltage. It is shown in [8] that after approximately five pumping cycles, the flow rate of the TWM micropump (without obstacle) finds a harmonic variation, which is in accordance with the harmonic deformation of the microchannel wall. The current study indicates that the same matter also holds for the proposed micropump (with obstacle). In other words, after approximately 5 cycles, the flow rate is the same for every cycle in the case of the TWM micropump, and for every other cycle in the case of the proposed micropump. Consequently, using the results of four cycles (which we selected to be the 7th to 10th cycles) seems to be adequate for comparing the flow rate of the two micropumps. The variation of the flow rate at the outlet of with- and without-obstacle micropumps is presented in Figure 6 for the 7th to 10th cycle of pumping. The first half of each cycle, where the flow rate curve is positive, corresponds to the pump phase in which fluid exits the microchannel from both the inlet and outlet, while in the second half of each cycle, the flow rate curve is negative, and fluid enters from both ends of the microchannel to accomplish the supply phase. The cyclic average of the outlet flow rate of each micropump is obtained by taking the average of its respective outlet flow rate (the curves
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
1265
Figure 7: The cyclic average of outlet flow rates of micropumps.
Figure 9: The fluid pressure contour for the with-obstacle micropump (in Pa).
Table 4: Cyclic average of the outlet flow rate of the 7th pumping cycle for different levels of element numbers. Number of elements 6 410 10 240 14 720 20 560 32 660
Figure 8: The fluid pressure contour for the without-obstacle micropump (in Pa).
of Figure 6) over one cycle. The averaging results for withoutand with-obstacle micropumps are presented in Figure 7 for seventh to tenth cycles. It can be seen from this figure that adding obstacles increases the cyclic average flow rate of the micropump from around 1 nL/s to around 6.5 nL/s, which is quite significant. This proves the efficiency of the proposed geometrical modification in increasing the net outlet flow rate of the Traveling Wave Micropump (TWM). The role of obstacles in the proposed micropump is to provide directional flow resistance (as explained in Section 2), which causes a better performance in terms of one-directional flow (from inlet to outlet). For example, in the 7th cycle, the total volume of fluid moved (i.e. back and forth) is 2.59 µL/s for the without-obstacle and 2.35 µL/s for the with- obstacle micropump. However, the cyclic average of the outlet flow rate, which is a measure of the total volume of fluid moved in the desired direction (from inlet to outlet), is 1.03 µL/s for the without-obstacle and 6.79 µL/s for the withobstacle micropump. It can be seen that in the without-obstacle
The outlet flow rate (nL/S) 5.48 6.07 6.62 6.79 6.79
micropump, a higher total volume of fluid is moved, but due to the presence of obstacles, the outlet flow rate of the with obstacle micropump is significantly higher. The FEM result for the cyclic average of the outlet flow rate of the 7th pumping cycle is presented in Table 4, for five different levels of element number. As this table indicates, discretizing the problem domain by 20 560 elements is sufficient for convergence of the simulations. A similar examination in time showed that a time step size of 1t = 1.25 × 10−7 (80 time steps per period of microchannel wall oscillation) suffices for obtaining a convergent result. The three-dimensional contour plots of the microchannel fluid pressure at different instances of the tenth cycle, for without- and with-obstacle micropumps, are presented in Figures 8 and 9, respectively. The maximum pressure difference generated in the microchannels of without- and with-obstacle micropumps is compared in Figure 10, for different instances of the tenth cycle. According to this figure, the maximum pressure difference is greater for with-obstacle micropump at all time instances. 5. Conclusion In this paper, a simple geometric modification to the traveling wave micropump has been proposed by placing two asymmetric trapezoid obstacles along the centerline of the micropump inlet and outlet. This modification is inspired by the action mechanism of the obstacle-type valveless micropump. The effect of this modification is examined through a verified finite element simulation of the proposed micropump. Fluid–structure interaction is considered in the simulation, and the variational multiscale method is used for fluid formulation in the ALE description. It is shown that adding
1266
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
Figure 10: The maximum pressure difference generated in the microchannels.
obstacles leads to approximately seven times increase in the outlet flow rate of the traveling wave micropump. References [1] Stemme, E. and Stemme, G. ‘‘A valveless diffuser/nozzle-based fluid pump’’, Sensors Actuators A, 39, pp. 159–167 (1993). [2] Koch, M., Harris, N., Evans, A.G.R., White, N.M. and Brunnschweiler, A. ‘‘A novel micromachined pump based on thick-film piezoelectric actuation’’, Sensors Actuators A, 70, pp. 98–103 (1998). [3] Ogawa, J., Kanno, I., Kotera, H., Wasa, K. and Suzuki, T. ‘‘Development of liquid pumping devices using vibrating microchannel walls’’, Sensors Actuators A, 152, pp. 211–218 (2009). [4] Lee, C.J., Sheen, H.J., Tu, Z.K., Lei, U. and Yang, C.Y. ‘‘A study of PZT valveless micropump with asymmetric obstacles’’, Microsyst. Technol., 15, pp. 993–1000 (2009). [5] Sheen, H.J., Hsu, C.J., Wu, T.H., Chang, C.C., Chu, H.C., Yang, C.Y. and Lei, U. ‘‘Unsteady flow behaviors in an obstacle-type valveless micropump by micro-PIV’’, Microfluid. Nanofluid., 4, pp. 331–342 (2008). [6] Pan, L.S., Ng, T.Y., Liu, G.R., Lam, K.Y. and Jiang, T.Y. ‘‘Analytical solutions for the dynamic analysis of a valveless micropump—a fluid-membrane coupling study’’, Sensors Actuators A, 93, pp. 173–181 (2001). [7] Zhou, Y. and Amirouche, F. ‘‘Study of fluid damping effects on resonant frequency of an electromagnetically actuated valveless micropump’’, Int. J. Adv. Manuf. Technol., 45, pp. 1187–1196 (2009).
[8] Afrasiab, H., Movahhedy, M.R. and Assempour, A. ‘‘Fluid–structure interaction analysis in microfluidic devices: a dimensionless finite element approach’’, Int. J. Numer. Methods Fluids, doi:10.1002/fld.2592 (2011). [9] Afrasiab, H., Movahhedy, M.R. and Assempour, A. ‘‘Finite element and analytical fluid–structure interaction analysis of the pneumatically actuated diaphragm microvalves’’, Acta Mech., doi:10.1007/s00707-0110508-9 (2011). [10] Donea, J. and Huerta, A., Finite Element Methods for Flow Problems, John Wiley & Sons, Chichester (2003). [11] Hsu, M.C., Bazilevs, Y., Calo, V.M., Tezduyar, T.E. and Hughes, T.J.R. ‘‘Improving stability of stabilized and multiscale formulations in flow simulations at small time steps’’, Comput. Methods Appl. Mech. Engrg., 199, pp. 828–840 (2010). [12] Stein, K., Tezduyar, T.E. and Benney, R. ‘‘Mesh moving techniques for fluid–structure interactions with large displacements’’, J. Appl. Mech., 70, pp. 59–63 (2003). [13] Ha, D.H., Phan, V.P., Goo, N.S. and Han, C.H. ‘‘Three-dimensional electrofluid–structural interaction simulation for pumping performance evaluation of a valveless micropump’’, Smart Mater. Struct., 18, pp. 104–111 (2009). [14] Fan, B., Song, G. and Hussain, F. ‘‘Simulation of a piezoelectrically actuated valveless micropump’’, Smart Mater. Struct., 14, pp. 400–405 (2005).
Hamed Afrasiab received his B.S. and M.S. degrees in Mechanical Engineering from Sharif University of Technology, Iran in 2004 and 2006, respectively, and is now a Doctoral degree student of the same subject at the same university. His main research area includes the finite element modeling of fluid–structure interaction in micromechanical systems. Mohammad Reza Movahhedy received his B.S. degree from University of Tehran, Iran in 1988, his M.S. degree from the University of Waterloo, Canada in 1994, and his Ph.D. degree from the University of British Columbia, Canada in 2000, all in Mechanical Engineering. He is currently Professor in the Department of Mechanical Engineering at Sharif University of Technology, Iran. His research interests are FEM simulation of metal cutting/forming processes, machine tool dynamics, mechanics of machining processes, experimental modal analysis and computer aided tolerancing. Ahmad Assempour received his B.S. and M.S. degrees from Tehran Polytechnic, Iran in 1979 and 1985, respectively, and his Ph.D. degree from Oklahoma State University, USA in 1989, all in Mechanical Engineering. He is currently Professor in the Department of Mechanical Engineering at Sharif University of Technology, Iran. His research interests are plasticity and mechanics of metal forming, modeling of sheet and bulk forming and die design in auto stamping industries.
Sharif University of Technology Scientia Iranica Transactions B: Mechanical Engineering www.sciencedirect.com
Proposal of a new design for valveless micropumps H. Afrasiab, M.R. Movahhedy ∗ , A. Assempour Center of Excellence in Design, Robotics and Automation, Department of Mechanical Engineering, Sharif University of Technology, Tehran, P.O. Box 11155-9567, Iran Received 7 February 2011; revised 9 August 2011; accepted 8 October 2011
KEYWORDS Valveless micropump; Finite element simulation; Fluid–Structure Interaction (FSI) analysis; Arbitrary Lagrangian–Eulerian (ALE) method.
Abstract A new design for a valveless micropumping device has been proposed that integrates two existing pumping technologies, namely, the wall induced traveling wave and the obstacle-type valveless micropump. The liquid in the microchannel is transported by generating a traveling wave on the channel, while the placing of two asymmetric trapezoid obstacles, along the centerline of the channel inlet and outlet, leads to a significant (up to seven times) increase of the net flow rate of the device. The effectiveness of this innovative design has been proved through a verified three-dimensional finite element model. Fluid–Structure Interaction (FSI) analysis is performed in the framework of an Arbitrary Lagrangian–Eulerian (ALE) method. © 2012 Sharif University of Technology. Production and hosting by Elsevier B.V. All rights reserved.
1. Introduction The growing requirement of controlled fluid transfer in many MEMS and BioMEMS applications has made the development of microfluidic systems an area of increasing interest in recent years. Dispensing drugs and other therapeutic agents, fuel delivery in miniaturized fuel cells and provision of refrigeration liquid in the cooling of microelectronic systems are some examples of these applications. Micropumps are essential components of microfluidic systems. Several different types of micropump have been developed up to now, valveless piezoelectric micropumps of which are in wide practical use amongst which valveless piezoelectric micropumps are in wide practical use, due to their simple structures, high flow rate output and ability to conduct particles such as red blood cells, polymers and proteins. These kinds of micropump are driven by a piezoelectric element bounded to a flexible membrane, and are equipped with a
∗
Corresponding author. E-mail address: [email protected] (M.R. Movahhedy).
1026-3098 © 2012 Sharif University of Technology. Production and hosting by Elsevier B.V. All rights reserved. Peer review under responsibility of Sharif University of Technology. doi:10.1016/j.scient.2011.11.023
pressure chamber for generating a pressure difference between the inlet and outlet, and no interior moving mechanical part is present in them. Instead, flow directing components like nozzle/diffuser elements or asymmetric obstacles are utilized to control rectification and obtain the net volume of the fluid flow [1,2]. Though piezoelectric micropumps can transport liquid at a high flow rate and with excellent controllability, the presence of the pressure chamber makes further miniaturization difficult and sometimes impractical [3]. Seeking a remedy for this problem, recently a mechanical micropumping device was fabricated and presented in [3], which uses a vibrating microchannel wall to drive the fluid flow while no pressure chamber is present. As reported in [3], on-generating a traveling wave on the channel wall, clear liquid flow could be induced within the microchannel. The transportation of the liquid by a traveling wave is controlled by changing the magnitude and frequency of the generated wave. A range of frequencies between 50 and 600 kHz were considered in [3] for oscillation of the microchannel wall for a maximum oscillation magnitude of 142.6 nm. A net outlet flow rate of 1.20 nL/s for the vibration frequency of 100 kHz was reported in [3] for a microchannel of 200 µm width and 300 µm height. The simple structure of this Traveling Wave Micropump (TWM) makes it advantageous for many microfluidic systems. In this paper, we propose a simple geometric modification to the traveling wave micropump so as to increase its net outlet flow rate. This modification is inspired by the action mechanism of the obstacle-type valveless micropump. The paper is organized as follows. The working principle of the new micropumping device is described in Section 2. The finite element modeling procedure is presented in Section 3. In Section 4, the results of FEM simulations are presented
1262
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
Figure 1: Working principle of the obstacle-type valveless micropump [5].
micropump. Now, we place two trapezoid obstacles along the centerline of the microchannel, one in the channel inlet downstream of the electrodes array and the other in the channel outlet upstream of the electrodes array, as shown in Figure 2. These obstacles form two equivalent planar nozzle/diffuser rectifying elements, similar to those described for the obstacletype valveless micropumps above. Consequently, a directiondependent flow resistance is provided, and an increase in the net flow rate of the microchannel is expected based on the previous experience with obstacle-type micropumps. In the next sections, we examine the effect of the proposed geometric modification on the net outlet flow rate of the microchannel, using a finite element simulation of both with- and withoutobstacle microchannels. 3. Finite element method simulation
Figure 2: A schematic of the proposed micropump.
and discussed. Finally, the concluding remarks are reviewed in Section 5. 2. Working principle of the new micropump In order to shed light on the action mechanism of the proposed micropump, first we briefly review the working principle of the obstacle-type valveless micropump presented e.g. in [4,5]. This kind of micropump is composed of a pump chamber with an oscillating diaphragm and two nozzle/diffuser flow rectifying elements, one directed from the inlet to the pump chamber, and the other from the pump chamber to the outlet, as illustrated in Figure 1. Each element is formed by two equivalent regions between the symmetrical trapezoid obstacle and the side-walls of the microchannel. Actuating the diaphragm causes a continual periodic increase and decrease in the volume of the pumping chamber (Figure 1). When the volume of the pumping chamber increases, the pressure in the chamber decreases, and more fluid enters through the nozzle-diffuser element on the left (inlet) relative to that on the right (outlet). This is because the element on the left acts as a diffuser, which poses less flow resistance than the nozzle on the right. Conversely, when the volume of the pumping chamber decreases, more fluid exits through the element on the right, which now acts as a diffuser. This results in a net pumping action from left to right, in Figure 1, as the diaphragm vibrates up and down. Now, we turn our attention to the new micropump. A schematic of this device is illustrated in Figure 2. In this device, the fluid flow is driven by the vibrating microchannel wall and actuated by a piezoelectric layer. The top ceiling wall of the microchannel is composed of a silicon wafer whose surface is covered by a piezoelectric PZT thin film. Nine separate top electrodes are attached to the PZT film. The traveling wave is induced on the ceiling wall by applying sinusoidal voltages, with a phase difference of 2π /3, to the adjacent top electrodes. When a channel wall oscillates in the form of traveling waves, peristaltic motion is induced in the liquid beneath the microchannel wall [3]. After a period of oscillation, the fluid moves slightly forward from the initial position, due to its viscosity, and the net liquid can be transported by repeating the period of motion. Up to this point, the pumping principle of the device is exactly similar to that of the traveling wave
In order to investigate the effect of obstacle addition on the net flow rate of the micropump, numerical simulations of both with- and without-obstacle microchannels were carried out in this section, using the finite element method. As stated earlier, the fluid flow is driven by the vibrating microchannel wall. Since the fluid also plays a role in resisting the wall vibration, the solid wall vibration and fluid flow are coupled, and a full Fluid–Structure Interaction (FSI) analysis is required in order to investigate the performance characteristics of these devices [6,7]. For this purpose, the FSI analysis is performed, using an Arbitrary Lagrangian–Eulerian (ALE) approach in the framework of the finite element method. The simulations are performed by means of a verified in-house developed FSIALE code, which has been used previously, e.g. in [8,9]. The finite element modeling details are described in the following subsections. 3.1. The fluid flow governing equations The law of conservation of momentum and the continuity equation for a viscous Newtonian incompressible flow are Navier–Stokes equations. In the ALE approach to fluid–structure interaction problems, these equations are written in the ALE description to account for deformation of the fluid mesh in the FSI interface. These equations, along with appropriate initial and boundary conditions, are seen in the ALE description (see e.g. [10])
ρf
∂ v − 2µf ∇.∇ s v + ρ v − vm .∇ v + ∇ p ∂t χ = ρ bf
in Ω f × (0, T ) ,
(1a)
∇.v = 0 in Ω f × (0, T ) ,
(1b)
v (x, 0) = v0
(1c)
on Ω × {0} , f
v (x, t ) = vD (x, t )
on
f ΓD
× (0, T ) ,
n .σ = −pn + 2µf n .∇ v = t f
f
f
f
s
f
on
(1d) f ΓN
× (0, T ) ,
(1e)
where ρf , µf , v and p represent the fluid density, viscosity, velocity and pressure, respectively. vm is the mesh velocity and bf denotes the body force vector. The fluid domain is specified by Ω f , while the time interval of interest is between (0, T ). The ∂(·)
time derivative in the ALE frame is represented by ∂ t . The χ
symmetric tensor, ∇ s v = ∇ v + (∇ v)T /2, is called the rate of deformation (or strain rate) tensor. nf is the exterior normal to the fluid boundary, Γ f , and v0 is the initial condition for
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
1263
f
the velocity field. vD denotes the prescribed velocity on the ΓD portion of the boundary, and tf is the boundary traction on the f complementary portion, ΓN . Since the working frequency of the solid wall is very high (between 50 and 600 kHz), the time step size required for the finite element analysis of these micropumps is very small. Thus, in order to avoid the so-called small time step instability in the finite element analysis of the fluid part, we used a residualbased variational multiscale formulation for incompressible Navier–Stokes equations, in combination with stabilization parameters, as recommended in [11], for weak formulation of the fluid flow. 3.2. The solid governing equations In a standard Lagrangian description, the law of conservation of linear momentum for a solid continuum may be expressed with respect to spatial coordinate, x, as:
∂ 2u = ∇.σs + ρs bs in Ω s × (0, T ) , ∂t2 u (x, 0) = u0 (x) in Ω s × {0} , ∂u (x, 0) = u˙ 0 (x) in Ω s × {0} , ∂t u (x, t ) = uD (x, t ) in ΓDs × (0, T ) , ρs
(2a)
Figure 3: The sequential fluid–structure coupling algorithm.
ρf vm .∇ v takes the convection of mesh moving points into
(2d)
account. The solid structure is formulated and solved in the Lagrangian description, as indicated earlier in this paper. Different techniques of ALE mesh moving are presented in the literature. Here, we use an approach that is based on solving the linear elasticity equations [12]. In this approach, the equation governing the displacement of fluid mesh nodes can be written as:
σ .n = t on × (0, T ) , (2e) where ρs is the solid density, and u represents its displacement
∇.σ = 0 in Ω f × (0, T ) , (5) where σ is the Cauchy stress tensor. For each boundary, i, a
s
s
s
(2b) (2c)
ΓNs
field, whereas the body forces are given by vector bs . The symmetric second order tensor, σs , denotes the Cauchy stress tensor, and the solid domain is represented by Ω s . uD denotes the prescribed velocity on the ΓDs portion of the boundary, and ts is the boundary traction on the complementary portion, ΓNs . The outer unit normal to the solid boundary is specified by vector ns . Since the solid membrane deformation is small in most microfluid manipulating devices [8], we assume that the structural part of the FSI problem is governed by linear elasticity. For a linear elastic solid, we have:
σs = λs (∇.u) I + 2µs ε (u) , (3) T in which ε (u) = ∇ u + (∇ u) /2 is the strain tensor, and µs and λs are the Lamé coefficients. 3.3. Fluid–structure interface conditions Fluid–solid interface conditions consist of kinematic and dynamic constraints specified as follows on the FSI interface, ∂ Ω fs :
∂ us = vf ∂t
on ∂ Ω fs × (0, T ) ,
σs .n + σf .n = 0 on ∂ Ω fs × (0, T ) ,
(4a) (4b)
with n being the outer normal at the solid boundary. 3.4. The mesh motion algorithm In the Arbitrary Lagrangian–Eulerian (ALE) approach to fluid–structure interaction problems, two key ingredients are needed. The first is a technique to move the fluid mesh, so that it can track structure motion at the FSI interface. This deformation is extended into the fluid interior field, in order to avoid excessive distortion of the fluid mesh at the interface. The second ingredient is a flow formulation written in the ALE framework, as presented in Eq. (1). The time derivative in Eq. (1a) is taken within the ALE description and the term,
Dirichlet boundary condition: ui = ubi ,
(6) b
may be given. u is the boundary displacement vector that is either given a priori or computed by the solid structure equations. 3.5. Fluid–structure coupling method In this work, a partitioned strong coupling approach is employed that uses a sequential iterative algorithm, as illustrated in Figure 3. This iterative scheme for the coupled system of fluid and structure is repeated at each time step until convergence is reached, i.e. after the norms of relative change in the field variables between two consecutive iterations are smaller than the specified convergence tolerance. 3.6. Boundary and initial conditions The boundary conditions are as follows: Along the lateral and bottom walls of the microchannel and obstacles (fixed boundaries), the usual ‘‘no slip’’ boundary condition is prescribed for the fluid. Along the top wall of the microchannel, which is also the bottom of the microchannel ceiling (FSI boundary), the fluid–structure interface conditions given in Sections 3 and 4 are applied. The motion of the microchannel ceiling is prescribed on its upper face. The initial condition is zero flow velocity and structure (microchannel ceiling) displacement. 3.7. Verification of the numerical approach In order to verify the accuracy of the finite element approach, a portion of the without-obstacle microchannel (containing all electrodes, as well as 1 mm downstream and 1 mm upstream of the electrodes array) was simulated in three dimensions, for 20 complete periods of oscillation and for different values of oscillation frequency. 14 040 and 3120 tri-quadratic elements
1264
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266 Table 1: Without-obstacle microchannel geometry. Silicon wafer thickness (µm) Piezoelectric film thickness (µm) Distance between single electrodes centers (µm) Excitation voltage (V) Maximum oscillation amplitude of electrodes for frequency of 100 kHz (nm) Microchannel width (µm) Microchannel height (µm) Microchannel length (mm)
Table 3: With-obstacle microchannel geometry. 6.0 2.5 400 20 142.6
Distance between two obstacles (L) (mm) Obstacles’ length (L1 ) (µm) Microchannel width (W ) (µm) Obstacles’ width (W1 ) (µm) Obstacles’ divergence angle (α )
3.2 500 200 25 5°
200 300 5.2
Table 2: Solid physical properties. Material Silicon Piezoelectric
Density (kg/m3 ) 2300 7800
Elastic modulus (GPa) 170 64
Poisson’s ratio
Figure 5: A schematic top view of the micropump geometry.
0.215 0.31
Figure 4: Frequency dependence of the mean flow velocity.
were used to discretize the fluid and solid parts, respectively. We used 80 time steps per period of oscillation in all these simulations. The without-obstacle microchannel geometry is presented in Table 1. The working fluid is water with density of 1000 kg/m3 and viscosity of 0.001 kg/(ms). The solid physical properties are listed in Table 2. In each case, the mean flow velocity was calculated in the microchannel outlet. Figure 4 compares the results of these simulations with the experimental data reported in [3]. It can be seen that the finite element results are in good agreement with experimental data, which demonstrates the capability of the employed numerical model. 4. Results and discussion Once the numerical model is verified, it can be used for studying the performance characteristics of the proposed micropump. To this end, in this section, a simulation similar to that of Section 3 was carried out for the with-obstacle micropump, for wall oscillation frequency of 100 kHz. 16 448 and 4112 tri-quadratic elements were used to discretize the fluid and solid parts, respectively. A schematic top view of the micropump geometry is indicated in Figure 5, and detailed dimensions of the obstacles and microchannel are listed in Table 3. The dimensions of the obstacles are adopted from those presented in [4]. Other geometric dimensions of the withobstacle micropump that are not given in this table are the same as those of the without-obstacle micropump presented in Table 1. It is to be noted that due to the presence of the
Figure 6: Outlet flow rate variation for without- and with-obstacle micropumps.
obstacles, more energy is needed for suction and discharge of fluid into and from the micropump. Consequently, in the new micropump, a higher voltage is needed for deforming the microchannel wall that controls the fluid flow. An electrofluid–structural interaction simulation (see e.g. [13,14]) is needed for determination of the required increase in the operating voltage. It is shown in [8] that after approximately five pumping cycles, the flow rate of the TWM micropump (without obstacle) finds a harmonic variation, which is in accordance with the harmonic deformation of the microchannel wall. The current study indicates that the same matter also holds for the proposed micropump (with obstacle). In other words, after approximately 5 cycles, the flow rate is the same for every cycle in the case of the TWM micropump, and for every other cycle in the case of the proposed micropump. Consequently, using the results of four cycles (which we selected to be the 7th to 10th cycles) seems to be adequate for comparing the flow rate of the two micropumps. The variation of the flow rate at the outlet of with- and without-obstacle micropumps is presented in Figure 6 for the 7th to 10th cycle of pumping. The first half of each cycle, where the flow rate curve is positive, corresponds to the pump phase in which fluid exits the microchannel from both the inlet and outlet, while in the second half of each cycle, the flow rate curve is negative, and fluid enters from both ends of the microchannel to accomplish the supply phase. The cyclic average of the outlet flow rate of each micropump is obtained by taking the average of its respective outlet flow rate (the curves
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
1265
Figure 7: The cyclic average of outlet flow rates of micropumps.
Figure 9: The fluid pressure contour for the with-obstacle micropump (in Pa).
Table 4: Cyclic average of the outlet flow rate of the 7th pumping cycle for different levels of element numbers. Number of elements 6 410 10 240 14 720 20 560 32 660
Figure 8: The fluid pressure contour for the without-obstacle micropump (in Pa).
of Figure 6) over one cycle. The averaging results for withoutand with-obstacle micropumps are presented in Figure 7 for seventh to tenth cycles. It can be seen from this figure that adding obstacles increases the cyclic average flow rate of the micropump from around 1 nL/s to around 6.5 nL/s, which is quite significant. This proves the efficiency of the proposed geometrical modification in increasing the net outlet flow rate of the Traveling Wave Micropump (TWM). The role of obstacles in the proposed micropump is to provide directional flow resistance (as explained in Section 2), which causes a better performance in terms of one-directional flow (from inlet to outlet). For example, in the 7th cycle, the total volume of fluid moved (i.e. back and forth) is 2.59 µL/s for the without-obstacle and 2.35 µL/s for the with- obstacle micropump. However, the cyclic average of the outlet flow rate, which is a measure of the total volume of fluid moved in the desired direction (from inlet to outlet), is 1.03 µL/s for the without-obstacle and 6.79 µL/s for the withobstacle micropump. It can be seen that in the without-obstacle
The outlet flow rate (nL/S) 5.48 6.07 6.62 6.79 6.79
micropump, a higher total volume of fluid is moved, but due to the presence of obstacles, the outlet flow rate of the with obstacle micropump is significantly higher. The FEM result for the cyclic average of the outlet flow rate of the 7th pumping cycle is presented in Table 4, for five different levels of element number. As this table indicates, discretizing the problem domain by 20 560 elements is sufficient for convergence of the simulations. A similar examination in time showed that a time step size of 1t = 1.25 × 10−7 (80 time steps per period of microchannel wall oscillation) suffices for obtaining a convergent result. The three-dimensional contour plots of the microchannel fluid pressure at different instances of the tenth cycle, for without- and with-obstacle micropumps, are presented in Figures 8 and 9, respectively. The maximum pressure difference generated in the microchannels of without- and with-obstacle micropumps is compared in Figure 10, for different instances of the tenth cycle. According to this figure, the maximum pressure difference is greater for with-obstacle micropump at all time instances. 5. Conclusion In this paper, a simple geometric modification to the traveling wave micropump has been proposed by placing two asymmetric trapezoid obstacles along the centerline of the micropump inlet and outlet. This modification is inspired by the action mechanism of the obstacle-type valveless micropump. The effect of this modification is examined through a verified finite element simulation of the proposed micropump. Fluid–structure interaction is considered in the simulation, and the variational multiscale method is used for fluid formulation in the ALE description. It is shown that adding
1266
H. Afrasiab et al. / Scientia Iranica, Transactions B: Mechanical Engineering 18 (2011) 1261–1266
Figure 10: The maximum pressure difference generated in the microchannels.
obstacles leads to approximately seven times increase in the outlet flow rate of the traveling wave micropump. References [1] Stemme, E. and Stemme, G. ‘‘A valveless diffuser/nozzle-based fluid pump’’, Sensors Actuators A, 39, pp. 159–167 (1993). [2] Koch, M., Harris, N., Evans, A.G.R., White, N.M. and Brunnschweiler, A. ‘‘A novel micromachined pump based on thick-film piezoelectric actuation’’, Sensors Actuators A, 70, pp. 98–103 (1998). [3] Ogawa, J., Kanno, I., Kotera, H., Wasa, K. and Suzuki, T. ‘‘Development of liquid pumping devices using vibrating microchannel walls’’, Sensors Actuators A, 152, pp. 211–218 (2009). [4] Lee, C.J., Sheen, H.J., Tu, Z.K., Lei, U. and Yang, C.Y. ‘‘A study of PZT valveless micropump with asymmetric obstacles’’, Microsyst. Technol., 15, pp. 993–1000 (2009). [5] Sheen, H.J., Hsu, C.J., Wu, T.H., Chang, C.C., Chu, H.C., Yang, C.Y. and Lei, U. ‘‘Unsteady flow behaviors in an obstacle-type valveless micropump by micro-PIV’’, Microfluid. Nanofluid., 4, pp. 331–342 (2008). [6] Pan, L.S., Ng, T.Y., Liu, G.R., Lam, K.Y. and Jiang, T.Y. ‘‘Analytical solutions for the dynamic analysis of a valveless micropump—a fluid-membrane coupling study’’, Sensors Actuators A, 93, pp. 173–181 (2001). [7] Zhou, Y. and Amirouche, F. ‘‘Study of fluid damping effects on resonant frequency of an electromagnetically actuated valveless micropump’’, Int. J. Adv. Manuf. Technol., 45, pp. 1187–1196 (2009).
[8] Afrasiab, H., Movahhedy, M.R. and Assempour, A. ‘‘Fluid–structure interaction analysis in microfluidic devices: a dimensionless finite element approach’’, Int. J. Numer. Methods Fluids, doi:10.1002/fld.2592 (2011). [9] Afrasiab, H., Movahhedy, M.R. and Assempour, A. ‘‘Finite element and analytical fluid–structure interaction analysis of the pneumatically actuated diaphragm microvalves’’, Acta Mech., doi:10.1007/s00707-0110508-9 (2011). [10] Donea, J. and Huerta, A., Finite Element Methods for Flow Problems, John Wiley & Sons, Chichester (2003). [11] Hsu, M.C., Bazilevs, Y., Calo, V.M., Tezduyar, T.E. and Hughes, T.J.R. ‘‘Improving stability of stabilized and multiscale formulations in flow simulations at small time steps’’, Comput. Methods Appl. Mech. Engrg., 199, pp. 828–840 (2010). [12] Stein, K., Tezduyar, T.E. and Benney, R. ‘‘Mesh moving techniques for fluid–structure interactions with large displacements’’, J. Appl. Mech., 70, pp. 59–63 (2003). [13] Ha, D.H., Phan, V.P., Goo, N.S. and Han, C.H. ‘‘Three-dimensional electrofluid–structural interaction simulation for pumping performance evaluation of a valveless micropump’’, Smart Mater. Struct., 18, pp. 104–111 (2009). [14] Fan, B., Song, G. and Hussain, F. ‘‘Simulation of a piezoelectrically actuated valveless micropump’’, Smart Mater. Struct., 14, pp. 400–405 (2005).
Hamed Afrasiab received his B.S. and M.S. degrees in Mechanical Engineering from Sharif University of Technology, Iran in 2004 and 2006, respectively, and is now a Doctoral degree student of the same subject at the same university. His main research area includes the finite element modeling of fluid–structure interaction in micromechanical systems. Mohammad Reza Movahhedy received his B.S. degree from University of Tehran, Iran in 1988, his M.S. degree from the University of Waterloo, Canada in 1994, and his Ph.D. degree from the University of British Columbia, Canada in 2000, all in Mechanical Engineering. He is currently Professor in the Department of Mechanical Engineering at Sharif University of Technology, Iran. His research interests are FEM simulation of metal cutting/forming processes, machine tool dynamics, mechanics of machining processes, experimental modal analysis and computer aided tolerancing. Ahmad Assempour received his B.S. and M.S. degrees from Tehran Polytechnic, Iran in 1979 and 1985, respectively, and his Ph.D. degree from Oklahoma State University, USA in 1989, all in Mechanical Engineering. He is currently Professor in the Department of Mechanical Engineering at Sharif University of Technology, Iran. His research interests are plasticity and mechanics of metal forming, modeling of sheet and bulk forming and die design in auto stamping industries.