A transient solution for droplet deformation under electric fields

A transient solution for droplet deformation under electric fields

Jia Zhang, Jeffery D. Zahn, and Hao Lin Department of Mechanical and Aerospace Engineering, Rutgers, The State University of New Jersey, Piscataway, NJ 08854, USA
Department of Biomedical Engineering, Rutgers, The State University of New Jersey, Piscataway, NJ 08854, USA
July 5, 2019
Abstract

A transient analysis to quantify droplet deformation under DC electric fields is presented. The full Taylor-Melcher leaky dielectric model is employed where the charge relaxation time is considered to be finite. The droplet is assumed to be spheroidal in shape for all times. The main result is an ODE governing the evolution of the droplet aspect ratio. The model is validated by extensively comparing predicted deformation with both previous theoretical and numerical studies, and with experimental data. Furthermore, the effects of parameters and stresses on deformation characteristics are systematically analyzed taking advantage of the explicit formulae on their contributions. The theoretical framework can be extended to study similar problems, e.g., vesicle electrodeformation and relaxation.

thanks: Email address for correspondence: hlin@jove.rutgers.edu

I Introduction

When a liquid droplet suspended in another immiscible fluid is subject to an applied electric field, it undergoes deformation due to the electrostatic stresses exerted on the interface. Extensive research on this phenomenon has been conducted to study the deformation due to its relevance in a variety of industrial applications, including electrohydrodynamic atomization Wu and Clark (2008), electrohydrodynamic emulsification Kanazawa et al. (2008), and ink-jet printing Basaran (2002), among others. Historically, the deformation dynamics is divided into two regimes: electrohydrostatics (EHS) and electrohydrodynamics (EHD). In the first, EHS deformation, the droplet is idealized as a perfect conductor immersed in a perfect insulating fluid; or both of the fluids are treated as perfect dielectrics with no free charge Allan and Mason (1962); Taylor (1964); Miksis (1981); Sherwood (1988); Ha and Yang (2000); Dubash and Mestel (2007). For this case, the electric field only induces a normal electrostatic stress, which is balanced by surface tension, and the final equilibrium shape is always prolate. At the steady state, the hydrodynamic flow is usually absent. In the second, EHD deformation, both fluids are considered to be leaky dielectrics Taylor (1966); Torza et al. (1971); Ajayi (1978); Sherwood (1988); Vizika and Saville (1992); Feng and Scott (1996); Baygents et al. (1998); Feng (1999); Hirata et al. (2000); Ha and Yang (2000); Bentenitis and Krause (2005); Sato et al. (2006); Lac and Homsy (2007). For this case, when an electric field is applied, free charges accumulate on the droplet surface which induces a tangential electrostatic stress in addition to the normal one. Driven by this force, the fluids inside and outside the droplet present toroidal circulations and a viscous stress is generated in response to balance the tangential electrostatic stress Taylor (1966). The droplet deforms into either a prolate or an oblate spheroid shape depending on the specific electrical properties of the fluids. With different electrical properties, the effects of the electrostatic and hydrodynamic stresses on droplet deformation are distinctive.

This work focuses on a solution method for problems of the second kind, namely, EHD deformation. This type of problem is more challenging to solve. In the literature, all theoretical solutions were obtained largely under two specific assumptions: (i) The deformations are small. The analysis is performed by assuming that the equilibrium shape of the droplet only slightly deviates from sphericity. Solutions using this assumption can be found in Taylor (1966); Torza et al. (1971); Ajayi (1978). (ii) For large deformations, the shape is assumed to be spheroidal during the entire deformation process. Results using this assumption are given in Bentenitis and Krause (2005). When compared with experimental data, predictions from the small-deformation theories always quantitatively underpredict the aspect ratio especially when the deformation is large. In contrast, the large-deformation theory has a better agreement both qualitatively and quantitatively. In all of the above, the theoretical analysis only leads to solutions in the steady state. The Taylor-Melcher leaky dielectric model Melcher and Taylor (1969); Saville (1997); Zhang et al. (2011) with the assumption of instantaneous charge relaxation has always been used. On the other hand, the theoretical analysis of transient droplet deformation seems to attract less attention. Only Dubash and Mestel Dubash and Mestel (2007) developed a transient deformation theory for an inviscid, conducting droplet. This analysis, which solves a EHS deformation problem, is not applicable to study EHD deformations. In general, to fully solve the transient EHD problem, numerical simulations have been employed Baygents et al. (1998); Feng (1999); Hirata et al. (2000); Lac and Homsy (2007).

In this work, we present a transient analysis of droplet deformation under direct-current (DC) electric fields. Following Bentenitis and Krause Bentenitis and Krause (2005), we assume the droplet remains spheroidal in shape. The full Taylor-Melcher leaky dielectric model is employed where the charge relaxation time is considered finite. In this framework, instantaneous charge relaxation is treated as a special limiting case. This generalization allows direct comparison with experimental data which were usually obtained in fluids with very low conductivities Ha and Yang (2000). The main result is an ordinary differential equation (ODE) governing the evolution of the droplet aspect ratio. The availability of this equation allows us to explicitly analyze the effects of parameters and stresses on the deformation characteristics. The model is validated by extensively comparing predicted deformation with both previous theoretical and numerical studies, and with experimental data.

Ii Theory

A schematic of the problem configuration is shown in Fig. 1(a). An uncharged, neutrally buoyant liquid droplet of radius is suspended in another fluid, and is subject to an applied electric field of strength . We assume that the fluids are immiscible leaky dielectrics with constant electrical and mechanical properties. and are the electrical conductivity, permittivity, and fluid viscosity, and the subscripts and denote internal and external, respectively. Under the influence of an applied electric field, free charges accumulate at the interface, which induces droplet deformation and EHD flows both inside and outside the droplet. Taylor Taylor (1966) predicted that droplets may deform into prolate or oblate shapes depending on the electrical properties of the fluids. In the following analysis, we focus on developing a solution for prolate deformations, whereas a solution for oblate deformations can be pursued in a similar manner (not presented here).

We assume that the droplet remains spheroidal in shape throughout the process. This approximation is consistent with experimental observations by Ha and Yang Ha and Yang (2000) and Bentenitis and Krause Bentenitis and Krause (2005). Following Tayor Taylor (1964), Bentenitis and Krause Bentenitis and Krause (2005), and Dubash and Mestel Dubash and Mestel (2007), the natural coordinate system to analyze this problem is the prolate spheroidal coordinate system, and a schematic is shown in Fig. 1(b). The geometry is assumed to be axisymmetric about the axis, which aligns with the direction of the applied electric field. The spheroidal coordinates are related to the cylindrical coordinates through the equations:

(1)
(2)

Here is chosen to be the semi-focal length of the spheroidal droplet, and and are the major and minor semi-axis, respectively. The contours for constant are spheroids, and . The contours for constant are hyperboloids, and . The surface of the prolate spheroid is conveniently given as

(3)

For the derivation below, we further assume that the volume of the droplet is conserved. We subsequently obtain

(4)

Therefore, the droplet geometry is completely characterized by a single parameter, , which evolves in time along with deformation. The critical idea of the current analysis is to express all variables, e.g., the electric potential and the stream function in terms of .

Figure 1: (a) A schematic of the problem configuration. (b) The prolate spheroidal coordinate system.

In what follows, we will solve the electrical problem first, followed by a solution of the hydrodynamic problem. An ODE for is obtained by applying both the stress matching and kinematic conditions.

ii.1 The electrical problem

The electric potentials inside and outside the droplet obey the Laplace equation according to the Ohmic law of current conservation with uniform electrical conductivity:

(5)

The matching conditions at the interface are

(6)
(7)

Here denotes a jump across an interface, and and are the unit tangential and normal interfacial vector, respectively. is the surface charge density. Note that in Eq. (7), we have included the displacement current, . This term is particularly important for fluids with very low conductivities (for example, those used in Ref. Ha and Yang (2000)) such that the interfacial charging time becomes comparable to the deformation time. However, we have ignored the effect of surface charge convection which is shown to be small by numerical simulations Feng (1999). Equation (7) can be rewritten in terms of the electric potentials,

(8)

Here is a metric coefficient of the prolate spheroidal coordinate system. The displacement current consists of two parts, represented by the first three terms on the LHS of the above equation. The first two terms result from a change in the droplet shape, and the third results from the charging process as if the shape remains unchanged. Here, we assume that the first term is negligible compared with the other two, and Eq. (8) can be further simplified to be

(9)

Indeed a consistency check a posteriori justifies the simplification. Far away from the droplet surface the electric field is uniform,

(10)

We also require that remains finite at . For the initial condition, we assume both the electric potential and the normal component of the displacement vector are continuous:

(11)

Solutions for the electric potentials have been obtained previously without including the displacement current Taylor (1966); Bentenitis and Krause (2005). With its inclusion the approach is similar and the results are

(12)
(13)

Here, is a 1st-degree Legendre polynomial of the second kind. is the dimensionless semi-focal length. The coefficients and are determined by the interfacial matching conditions (6) and (9) which gives

(14)
(15)
(16)

Here and are the permittivity ratio and the conductivity ratio, respectively. is an electrical charging time. is a characteristic flow timescale used below in the hydrodynamic problem, and is the coefficient of surface tension. In the above equations, a dimensionless time has been used. In general, Eq. (15) needs to be integrated together with an ODE for to obtain and . However, in the limit of instantaneous-charge-relaxation time, , and Eq. (15) can be simplified to be

(17)

This result is equivalent to a solution employing the simplified boundary condition in place of Eq. (7).

The normal and tangential electrostatic stresses are given by,

(18)

where and are the normal and tangential electric fields, respectively. is a metric coefficient of the prolate spheroidal coordinate system. These stresses can be evaluated with the solutions (12) and (13), and will be used in the stress matching conditions below.

ii.2 The hydrodynamic problem

In the regime of low-Reynolds-number flow, the governing equation for the hydrodynamic problem can be rewritten in terms of the stream function, , as

(19)

Here, the expression for the operator can be found in Dubash and Mestel Dubash and Mestel (2007) and Bentenitis and Krause Bentenitis and Krause (2005). The stream function is related to the velocity components as

(20)

is a metric coefficient of the prolate spheroidal coordinate system. At the interface, and represent the tangential and normal velocities, respectively, and they are required to be continuous

(21)

In addition, we prescribe a kinematic condition relating the interfacial displacement to the normal velocity,

(22)

The total force on the interface resulting from the electrical stress, the hydrodynamic stress, and the surface tension should be balanced at every point. However, this constraint is impossible to satisfy exactly within the framework of spheroidal deformation. Various authors developed reduced stress-balance conditions instead Taylor (1964); Sherwood (1988); Bentenitis and Krause (2005); Dubash and Mestel (2007). Here we follow the integrated formulae proposed by Sherwood Sherwood (1988) and Dubash and Mestel Dubash and Mestel (2007)

(23)
(24)

Equations (23) and (24) represent a global balance of the tangential and normal stresses, respectively derived from energy principles. Here denotes the hydrodynamic stress, and are the two principal radii of the curvature, and the integration is carried over the interface.

The general solution to (19) was proposed by Dassios et al. Dassios et al. (1994) using the method of semi-separation:

(25)

Here and are Gegenbauer functions of the first and second kind, respectively. and are linear combinations of and . The detailed expressions for , , , and are found in Ref. Dassios et al. (1994). Interested readers are referred to Ref. Dassios et al. (1994) for further details. After considering that the far field is quiescent, and that the velocities remain finite at , the stream functions can be simplified to be

(26)
(27)

where and are unknown coefficients satisfying the relations , . In general, these coefficients are inter-dependent, and the full solution can be obtained only with the entire infinite series. Here we seek a truncated solution as an approximation,

(28)
(29)

Indeed, gives a functional form in confirming with that in Eq. (22), which can be rewritten as

(30)

This agreement in part validates the spheroidal shape assumption: the shape represents the leading mode in the infinite series.

Equations (21-24) are combined to solve for the five unknown variables, namely, , , , , and . Specifically, Eqs. (21) and (22) are first used to eliminate the , , ,

(31)
(32)
(33)

where , and . Further considering Eq. (23), we can express in terms of ,

(34)

where is the viscosity ratio. The detailed expressions of are found in the Appendix. This expression is inserted into Eq. (24) to obtain the final result, an ODE governing the evolution of the ,

(35a)
(35b)
(35c)

The detailed expressions of , and are also found in the Appendix. The coefficients and are given by Eqs. (14) and (15), respectively. is the electric capillary number. In Eq. (35a), the three terms in the numerator on the RHS represent the contributions from the normal stress, the tangential stress, and the surface tension, respectively. At the equilibrium, the balance of the three forces determines the final shape. The leading coefficients and arise exclusively from the electrostatic stresses, and can be used to estimate their respective influence on deformation. In the limit of instantaneous relaxation, and by considering Eqs. (14) and (17), and can be simplified to be

(36)
(37)

For this case, the evolution of is governed by a single timescale, . Once is obtained by solving the Eqs. (15) and (35a), the aspect ratio is calculated by the formula

(38)

Iii Comparison with previous results

In this section, we compare our model prediction extensively with results from previous work. The comparisons with theoretical/numerical results and experimental data are respectively presented in Secs. III.1 and III.2.

iii.1 Comparison with previous theories and simulation

We first consider the equilibrium shape, and compare our results with those from Bentenitis and Krause Bentenitis and Krause (2005). For this case, the LHS of Eq. (35a) is simply set to zero, resulting in the so called discriminating equation,

(39)

Here and are given by Eq. (36). is solved as a root(s) of this equation from which the equilibrium aspect ratio, , can be obtained. Equation (39) shows that the equilibrium shape is only determined by the dimensionless parameters and . A comparison with the theoretical prediction by Bentenitis and Krause Bentenitis and Krause (2005) is shown in Fig. 2. Note that in this earlier work, the authors solved for the equilibrium shape directly without obtaining the transient solution. A good agreement is observed, although a different stress matching condition has been used by Bentenitis and Krause Bentenitis and Krause (2005) [see Eqs. (38) and (45) therein].

Figure 2: The equilibrium aspect ratio as a function of electric capillary number. The parameters are , , and .

We next compare with the results from Dubash and Mestel Dubash and Mestel (2007). In this work, the authors developed a theoretical model, also with the spheroidal shape assumption, to predict the transient deformation of a conducting, inviscid droplet immersed in a viscous, nonconductive solution. This special consideration leads to significant simplifications: both the electric and hydrodynamic fields are absent within the droplet. In addition, at the equilibrium state (if one is permitted), the hydrodynamic flow outside the droplet is also quiescent, giving rise to the phenomenon termed EHS.

In our generalized framework, the solution for this case is simply achieved by setting and in Eqs. (35a) and (36). Note that directly leads to instantaneous charge relaxation. The resulting comparisons are shown in Fig. 3 in which the aspect ratio () is plotted as a function of time for four different electric capillary numbers (). For the two lower values of , the current model has excellent agreement with both the theoretical and numerical predictions by Dubash and Mestel Dubash and Mestel (2007) [Fig. 3(a)]. For these values, final equilibria are achieved. As increases [Fig. 3(b)], the deformation becomes unstable and an equilibrium shape is no longer possible. The rapid expansion with a sharp slope at the later stage preludes droplet breakup. For these two cases, the theoretical models still agree with each other, whereas some discrepancies exist with respect to the numerical simulation, in particular for =0.206. However, this discrepancy is in general only noticeable when the number is above and very close to the critical threshold of breakup (0.2044 for the case studied), due to a slight underprediction of the rate of deformation by the theoretical models. A similar trend is observed when comparing with the numerical simulation by Hirata et al. Hirata et al. (2000) (not shown). Overall, our model can serve as a good approximation to the numerical model which is considered more accurate.

Figure 3: The deformation of a conducting droplet in a highly viscous medium. (a) and 0.204. (b) and 0.21. The dimensionless time is defined as , where .

iii.2 Comparison with experimental data

The main source of experimental data comes from Ha and Yang Ha and Yang (2000). We also begin with an examination of the final aspect ratio when an equilibrium shape can be achieved. Figure 4 shows the equilibrium aspect ratio of a castor oil droplet immersed in silicone oil from Ha and Yang Ha and Yang (2000), as well as predicted by various models. The current prediction is shown as a solid line, whereas the results from first-order Taylor (1966) and second-order Ajayi (1978) theories are shown as dot-dashed and dashed lines, respectively. Following Lac and Homsy Lac and Homsy (2007), we rescale to best match Ajayi’s second-order correction. This rescaling is equivalent to adjusting the surface tension from used by Ha and Yang Ha and Yang (2000) (which is a fitting parameter in that work) to . The latter value is close to the lower bound, , measured by Salipante and Vlahovska Salipante and Vlahovska (2010). In addition, we use according to the measurements by Torza et al. Torza et al. (1971), Vizika and Saville Vizika and Saville (1992), and Salipante and Vlahovska Salipante and Vlahovska (2010), which is slightly different from the value of used by Lac and Homsy Lac and Homsy (2007). The results show good agreement between the current model and the experimental data. Most importantly, our theory correctly predicts a critical (0.244) for droplet breakup. In contrast, the small deformation theories can not capture this critical phenomenon.

We have also compared our theoretical prediction with the experimental data from Bentenitis and Krause Bentenitis and Krause (2005), which measured the equilibrium aspect ratio of a DGEBA droplet immersed in a PDMS solution. Since our result is in good agreement with the theoretical prediction in the same work (see Fig. 2), which in turn agrees well with the data, the comparison is not shown here for brevity.

Figure 4: The prediction from current model is compared with the small deformation theories Taylor (1966); Ajayi (1978) and experimental data Ha and Yang (2000). The parameters are , , and .

Next, we will compare the transient solution from our model with data from Ha and Yang Ha and Yang (2000). In Fig. 5(a), the data is extracted from Fig. 3 in the latter work, which captures the deformation of a water droplet in silicone oil. The droplet is fitted with an ellipse at every instant, based on which the aspect ratio is calculated. A 10% fitting error is estimated, and is shown as error bars in Fig. 5(a) [the same approach is adopted to extract the data presented in Figs. 5(b) and 6]. The model prediction is calculated with Eqs. (35a) and (36), and with , , , , , and all directly taken from Ha and Yang Ha and Yang (2000). For medium permittivity, we use following the measurements by Torza et al. Torza et al. (1971), Vizika and Saville Vizika and Saville (1992), and Salipante and Vlahovska Salipante and Vlahovska (2010). For surface tension, we use , which is consistent with the values reported by Torza et al. Torza et al. (1971) and Vizika and Saville Vizika and Saville (1992). In this case, the model is able to predict the deformation process with good quantitative accuracy. In Fig. 5(b), a similar comparison is shown for a water-ethanol droplet in silicone oil. The data is based on Fig. 4 in Ref. Ha and Yang (2000). For our calculation, , , , , , , and . Because the droplet is doped with polyvinylpyrrolidone (a polymer solution), the surface tension is not directly available, and is used as a fitting parameter instead to generate the best agreement between theory and data. The resulting value is , 11% higher than that for water/silicone oil which is used in Fig. 5(a).

Figure 5: Comparison of transient droplet deformation. (a) A water droplet in silicone oil. The parameters are , , , , , , , and . (b) A water-ethanol droplet in silicone oil. The parameters are , , , , , , , and .
Figure 6: Droplet deformation in the limit of extremely low conductivities. The parameters are , , , , , , , and . The best agreement between the data and the theory is found for . For reference, the dotted line shows the calculation according to the instantaneous-charge-relaxation (ICR) model.

In contrast to the regime of instantaneous charge relaxation examined in Fig. 5, Fig. 6 represents droplet deformation in the finite-charging-time regime. The data is extracted from Fig. 7 in Ref. Ha and Yang (2000). For this case, the droplet is made of castor oil, and is immersed in silicone oil. The extremely low conductivities of these media lead to a charging time (seconds) comparable to the deformation time, and the full model, Eqs. (35a)-(35c), has to be used. For our calculation, , , , , , , , and . Note that the values for the surface tension and the conductivity ratio follow the measurements by Torza et al. Torza et al. (1971), Vizika and Saville Vizika and Saville (1992), and Salipante and Vlahovska Salipante and Vlahovska (2010) which are believed to be more accurate than the original values of and given by Ha and Yang Ha and Yang (2000). In addition, the actual conductivity of silicone oil varies from to in the literature Park et al. (2009); Salipante and Vlahovska (2010); Kim et al. (2011). In Fig. 6, we show the calculation with three representative values within this range, namely, , , and . The best agreement is found for . For comparison, the calculation according to the instantaneous-charge-relaxation model [Eqs. (35a) and (36)] is also shown, and is denoted by ICR. This simplified model clearly overpredicts deformation by a significant degree.

In general, our model agrees well with experimental data in both steady and transient states, and for a large parametric range. These comparisons provide a strong validation for our model.

Iv The effects of stresses on deformation

In this section, we demonstrate the utility of our theoretical results by analyzing in-depth the governing equation. For simplicity, we focus on the regime of instantaneous relaxation, where and are given by Eq. (36). A main contribution of the current work is that Eq. (35a) clearly separates the effects by different forces. In the numerator of the RHS, the three terms represent respectively the effects of the normal stresses (both electrical and hydrodynamic), the tangential stresses (both electrical and hydrodynamic), and the surface tension. Furthermore, all the functions in this equation are positive (, , , ), such that the signs of and completely determine whether the normal and tangential stresses would promote or suppress deformation. Due to the inverse relationship between and the aspect ratio, [see Eq. (38)], a positive or indicates a positive contribution. Evidently, surface tension always resists deformation. Because and depend exclusively on the electrical properties in a simple manner [see Eq. (36)], their influences can be conveniently analyzed using a phase diagram shown in Fig. 7. The dashed and dotted lines correspond to and , respectively. These lines separate the phase space into three regimes, where and denote the normal and tangential stresses, and the superscripts and denote a positive or negative contribution to deformation, respectively. In addition, the solid line is obtained by solving for the root of Taylor’s discriminating function Taylor (1966), which separates the prolate (denoted by ) and oblate (denoted by ) regimes [this line can be equivalently obtained by looking for the steady-state solution of from Eq. (39)].

Figure 7 can be used to shed light on the physical processes governing deformation. First, the line for separates the and regimes, which corroborates with the previous results Taylor (1966); Lac and Homsy (2007). On this dividing line, the velocity field becomes zero, so does the tangential electrical stress. In Ref. Lac and Homsy (2007), the viscosity ratio has opposite effects on deformation in the and regimes. This behavior is clearly explained by Eq. (35a). Second, there is a small region within the oblate regime, namely, the area between the solid and dashed lines where is positive. This suggests that the normal stress still tends to stretch the droplet along the direction of the applied field. However, because is negative, the tangential stresses overcome the normal stresses, and stretch the droplet into an oblate shape. This new insight is not available from previous analysis or simulations.

Third, in the prolate regime where is always positive, the sign of leads to different deformation behavior. Figure 8 shows the equilibrium aspect ratio as a function of for three specific cases. Note that the new variable

(40)

In this new definition, corresponds to , and corresponds to . For all three cases, and . For , . We observe hysteresis, and approaches 1 rapidly in the upper brunch. The cases of and correspond to and , respectively. In general, as decreases, the deformation becomes weaker for comparable values. Most interestingly, for (), converges to a value less than 1 in the limit of . This means that even for the very large applied electric field strength, a finite equilibrium aspect ratio can be achieved. We emphasize this scenario is only possible in the regime. For large values, corresponding to large , the resistive effect from surface tension is negligible, and the only way to obtain a finite equilibrium aspect ratio is therefore by balancing the normal and tangential stresses. Since is positive, has to be negative.

Figure 7: Phase diagram for droplet deformation. Here the dotted line is calculated by satisfying . The solid line is calculated by solving for the root of Taylor’s discriminating function. The dashed line represents . and denote prolate and oblate deformation, respectively. and denote the effect of normal and tangential stresses, respectively, and a or sign denotes facilitating or suppressing, respectively.
Figure 8: The behavior of equilibrium droplet deformation in different regimes. For , ; , ; , . As , an equilibrium shape is only possible in the regime. Other parameters are and .

V Conclusions

In conclusion, we have developed a transient analysis to quantify droplet deformation under DC electric fields. The full Taylor-Melcher leaky dielectric model is employed where the charge relaxation time is considered finite. In this framework, instantaneous charge relaxation is treated as a special limiting case. The droplet is assumed to be spheroidal in shape for all times. The main result is an ODE governing the evolution of the droplet aspect ratio. The model is validated by extensively comparing predicted deformation with both previous theoretical and numerical studies, and with experimental data. In particular, the experimental results by Ha and Yang Ha and Yang (2000), which were obtained with extremely low medium conductivities are well captured by the simulation with the finite-time charge-relaxation model. The model is used to analyze the effects of parameters and stresses on the deformation characteristics. The results demonstrate clearly that in different regimes according to the sign of , the stresses contribute qualitatively differently to deformation. Last but not least, this work lays the foundation for the study of a more complex problem, namely, vesicle electrodeformation and relaxation. This problem is the pursuit of our future work.

Acknowledgements.
JZ and HL acknowledge fund support from an NSF award CBET-0747886 with Dr William Schultz and Dr Henning Winter as contract monitors.

Appendix A Appendixes

The functions in Eq. (34) are given in the following expressions:

(41)
(42)
(43)
(44)
(45)

The functions and in Eq. (35a) are given in the following expressions:

(46)
(47)
(48)
(49)
(50)

where