Effect of input voltage frequency on the distribution of electrical stresses on the cell surface based on single-cell dielectrophoresis analysis

Electroporation is defined as cell membrane permeabilization under the application of electric fields. The mechanism of hydrophilic pore formation is not yet well understood. When cells are exposed to electric fields, electrical stresses act on their surfaces. These electrical stresses play a crucial role in cell membrane structural changes, which lead to cell permeabilization. These electrical stresses depend on the dielectric properties of the cell, buffer solution, and the applied electric field characteristics. In the current study, the effect of electric field frequency on the electrical stresses distribution on the cell surface and cell deformation is numerically and experimentally investigated. As previous studies were mostly focused on the effect of electric fields on a group of cells, the present study focused on the behavior of a single cell exposed to an electric field. To accomplish this, the effect of cells on electrostatic potential distribution and electric field must be considered. To do this, Fast immersed interface method (IIM) was used to discretize the governing quasi-electrostatic equations. Numerical results confirmed the accuracy of fast IIM in satisfying the internal electrical boundary conditions on the cell surface. Finally, experimental results showed the effect of applied electric field on cell deformation at different frequencies.

parameters such as permittivity, conductivity, and material parameters combined with electric field frequency are applied. These modeling tools can be used for accurate prediction of the electric field-induced cell shape deformations 18 .
Several studies have been carried out on electric field-induced cell membrane permeabilization. For example, Neu et al. 19 carried out numerical simulations of radial electromechanical force for toroidal and cylindrical pores in the presence of a uniform electric field and calculated the electrical energy of the nanopores as a function of pore size. Molecular dynamics simulation was performed by Tieleman et al. 20 to show the ability of pore deformations under appropriate electromechanical forces. Some experimental results suggest that the poration of vesicles is highly dependent on the electric field exposure time. Salipante et al. 21 showed that vesicles would rupture in a weak electric field when stressed for a long time, while they survive a strong electric field for short durations. The effect of uniform electric fields on transmembrane potential has been studied numerically by Priva and Gowrisree 22 . Differences between the inner and outer membrane potential have been shown by the comparison of numerical results and analytical values. Needham and Hochmuth 23 studied the effect of critical external electric field strength on mechanical stretching of the cell membrane. Electroporated cells showed considerable changes in geometrical and electrical properties, which was due to membrane breakdown of lipid vesicles under an external electric field. Numerical and experimental studies on the effect of such cellular modifications were carried out by Oblak et al. 24 and were applied as a basis of an electric field mediated cell separation technique. A novel chip had been designed to calculate the dynamic changes in membrane permeability through electroporation in human cancer cells by He et al. 25 In this paper, dynamics of pore resealing were also investigated and the resealing time constants were calculated for different pulse treatments. In another research, the effect of a cell's volume fraction, medium conductivity, membrane conductivity, critical transmembrane potential, and cell orientation have been investigated on the effective conductivity of a suspension of electroporated cells 26 . The anisotropic nature of electropermeabilized cells was also evaluated using the finite element method and an analytical approach. Zudans et al. 27 represented a numerical calculation method for estimating the degree of cell electropermeabilization in inhomogeneous electric fields. In the aforementioned studies, the effect of electric field characteristics on cell permeabilization has been studied using different experimental and numerical methods such as molecular dynamics. However, the effect of cell deformation and cell membrane structural changes caused by electrical stresses has not been discussed. Electrical stress distribution over the cell surface and its effect on membrane permeabilization in the presence of an electric field is an imperative process that requires fundamental investigations and should not be overlooked.
The exact mechanism of pore formation across the plasma membrane in the presence of an electric field is not yet well understood. Pipet aspiration is another method used for cell membrane permeabilization. In this method, mechanical tension is applied to cells. Thus it seems that cell deformation caused by mechanical stresses plays an essential role in pore formation on the cell membrane. As mentioned earlier, cells experience such stresses in the presence of an electric field. To better understand the electroporation phenomena, it is necessary to study the effect of electric fields on cells. To this end, the present work studied the behavior of a red blood cell under the application of a non-uniform AC electric field. Electric stresses induced by an electric field acts on the cell surface and cause cell deformation, which may play a role in cell membrane structural changes and its permeability.
Despite the considerable interest in the field of electroporation, to the best of our knowledge, no research group has studied the effect of the applied electric field frequency on the electric stress distribution over a cell's surface and its subsequent effects on membrane permeabilization analytically. In this paper, the immersed interface method (IIM) was employed to study the electrical stresses generated on a cell membrane exposed to a non-uniform AC electric field. Electrical stresses were calculated on the inner and outer boundaries of the cell membrane using the Maxwell stress tensor (MST). To validate the numerical algorithm, the resultant of electric tension applied on the cell membrane was compared with the dielectrophoretic (DEP) force, which was calculated analytically using effective dipole moment (EDM) approximation. Finally, an experimental analysis was conducted to present the effect of electric stresses on cell deformation at different frequencies.

Theory
We considered a cell with the conductivity of σ p , permittivity of ε p , immersed in a dielectric fluid with ohmic conductivity of σ f , and permittivity of ε f . A non-uniform AC voltage is applied on the enclosure. The geometry of electrodes is presented in Fig. 1.
The well-known governing quasi-electrostatic equation is seen in Eq. 1, where ε is the complex permittivity that is defined as Eq. 2 with ω π = f 2 and f being the applied AC voltage frequency.
In the above equation, ϕ is the complex electric potential that is expressed as Eq. 3.
The boundary conditions for the enclosure, shown in Fig. 1, are as follow: The internal boundary condition (jump condition) on the cell surface is as Eq. 5.
The effect of electric stress acting on a cell's surface that is being subjected to an electric field is calculated by using MST. By neglecting the magnetic effects, MST can be calculated for any dielectric particle with ohmic losses using Eq. 6.
In the above equation, the electric field can be any function of time. According to the literature, when a cell is exposed to an AC electric field, the time scale of its motion is much larger than the electric field frequency. Hence, the time average of MST is calculated, and the time-dependent part is ignored. The electric field can be expressed as Eq. 7.
M where ⁎ E is the complex conjugate of E . Maxwell stress tensor for a cell with an ohmic loss can be written as Eq. 8.
The first part of Eq. 8 is the time average MST, and the second part is its time-dependent part, which goes to zero under the integration over the electric field time period. Thus, the time average of electric stresses acting over the cell surface can be obtained as Eq. 9, in which The integration of electric stresses over the cell surface is equal to the net dielectrophoretic (DEP) force exerted on a cell in the presence of the non-uniform electric field. Using the effective dipole moment (EDM) approximation, the time average of DEP force can be calculated using Eq. 10. where R is the cell radius and ω K( ) is the Clausius-Mossotti factor, and they depend on the dielectric properties of the cell, suspension media, and the frequency of the applied electric field. Equation 11 defines the Clausius-Mossotti factor for a 2D cell with a complex permittivity of ε p , immersed in a fluid with the permittivity of ε f .
at the frequency f.
The effective dipole moment is valid only when the cell size is negligible in comparison with the enclosure size or electric field nonuniformity. In such cases, the effect of the cell on the electric potential distribution and the electric field can be neglected. Otherwise, EDM approximation will be inaccurate. The real part of the Clausius-Mossotti factor determines the behavior of the cell in the presence of an electric field.

Solving Method
To calculate MST on the cell surface, the governing quasi-electrostatic equation in the presence of a cell was first solved. In previous studies, the existence of the cell has been mostly ignored in the computational domain and the cell is replaced by equivalent multipoles. In the present study, fast (augmented) IIM was used to implement internal electric boundary conditions on the cell surface. IIM is a modified finite difference (or finite element) method developed by Z. Li for solving partial differential equations (PDEs) involving interfaces and irregular domains.
In the current study, the governing quasi-electrostatic equation was derived using fast IIM. Fast IIM is a modified immersed interface method for solving PDEs, where the coefficients have a constant value in each subdomain. By considering a 2D elliptic equation involving interface Γ, β x y ( , ) is assumed to be constant in each domain (Eq. 12).
Internal boundary conditions (jump conditions) on the interface Γ are as follow: The discrete form of Eq. 12-a can be written as: In the above equation, C ij is the correction term. If the coefficient β is constant in Eq. 12, or is constant in each of the domains, or there is only the singular force on the interface, then the finite difference coefficients of Eq. 14 are simplified to the standard 5-point finite difference scheme. However, a correction term still needs to be added, which is due to a jump in the value of coefficient β or a source distribution along the interface. Thus, the discrete form of Eq. 12 can be written as Eq. 15, where Δ h is the discrete Laplacian operator and C ij is the correction term.
Based on the standard 5-point finite difference stencil centered at i j ( , ), if all five grid points are on the same side of the interface, point i j ( , ) is called a regular grid point. Otherwise, it is an irregular grid point. The correction term at regular points is zero; grid points that are away from interfaces. At irregular points, the correction term is calculated in terms of the first and the second surface derivatives of ν and ω at the control point X Y ( , ) on the interface using Eq. 16, where χ is the interface curvature.   In the above equations, the values of γ k are determined according to the applied finite difference stencil and . ξ k and η k are local coordinates of the grid points' finite-difference stencil centered at i j ( , ) in the system that its origin is the control point X Y ( , ) on the interface in the normal and tangential directions.
The augmented variable can be considered as = . The discrete values of ω and g on the interface at each of the control points are , respectively and n b is the number of control points. Since the correction term is a linear combination of G { } k and W { } k , the matrix form of Eq. 15 can be written as follow: The discrete form of jump condition (Eq. 13-a) is obtained as Eq. 19, where U n is the discrete value of ∂ ∂ u n / and V is the discrete values of ν on the interface.  In Eq. 20, C is the correction term. By solving the following system of equations, the coefficients γ { } k can be determined. In the case of an undetermined system of equation, singular value decomposition (SVD) must be applied. By using a weighted least square interpolation, the discrete form of jump condition was obtained as Eq. 22. where C is expressed as: Finally, the jump condition was obtained as Eq. 24. By solving Eqs. 18 and 24, the distribution of u can be calculated.
To examine the accuracy of fast IIM, Eq. 12 was solved using IIM for β = www.nature.com/scientificreports www.nature.com/scientificreports/ The numerical and exact solutions of the problem are plotted in Fig. 3. In Fig. 4, the numerical and the exact solutions are compared at lines = x 0 and = y 0. These figures show that the solution algorithm is accurate, and the internal boundary conditions are well satisfied. To calculate the electric stresses acting over the cell surface in the presence of an AC electric field, the electric field discontinuity on the cell membrane must be calculated accurately. These results show that IIM is an appropriate numerical method, which helps to handle the field parameter discontinuity in the solution domain using a uniform mesh.
By separating the real and the imaginary parts of Eq. 1, the following equation is obtained: If the dielectric properties of the cell are isotropic, it can be concluded that: Also, the jump condition on the cell surface can be presented as follow: By considering the augmented variables, g 1 and g 2 , as follow: The discrete form of Eq. 26 can be obtained as Eq. 29, where G 1 and G 2 are discrete values of g 1 and g 2 .
The discrete jump condition of Eq. 27 was derived as Eq. 30, where F 3 and F 4 are zero in this problem.
Finally, Eq. 31 can be solved to obtain the electric potential distribution in the enclosure in the presence of a cell.     www.nature.com/scientificreports www.nature.com/scientificreports/

Results and Discussion
Numerical solution. Demonstrating convergence in numerical analysis is the first step in solving problems.
The mesh convergence study, for different mesh sizes, is shown in Fig. 5. This figure displays the electrostatic potential at = x 25 μm (middle vertical line in Fig. 1). In addition to the mesh size, the numbers of interpolation and control nodes are also important for the convergence of IIM. The number of control nodes was determined according to the mesh size. Sixteen nodes were used for interpolation. According to Fig. 5, the mesh size = .
h 0 6 μm was shown appropriate for a good convergence. As illustrated in these figures, the electric potential on the boundary of the cell is continuous while the electric field, which is the electric potential gradient, is discontinuous. Figure 6 presents the electrostatic potential distribution and the electric field in the channel without the presence of a cell.
The value and the direction of the DEP force depend on the value and the sign of the real part of the Clausius-Mossotti coefficient. The dielectrophoretic spectrum of a cell, which shows the effect of the frequency on the real part of the coefficient, provides information about the cell's electrical structure. In addition to the frequency, the real part of the Clausius-Mossotti coefficient depends on the dielectric properties of the cell and the suspension buffer. In Fig. 7, the effect of frequency on the real part of the Clausius-Mossotti is shown.     www.nature.com/scientificreports www.nature.com/scientificreports/ As can be seen in Fig. 7, the behavior of the cell changes in the frequency range from 1 kHz to 10 MHz. Thus, to determine the effect of frequency on cell behavior, the equations were solved for this range of frequencies. The electric potential and the electric field distribution in the presence of the cell at the frequencies of 10 kHz and 10 MHz were calculated. Different effects of a cell on the distribution of the electric potential and the electric field at the two discussed frequencies are presented in Figs. 8 and 9.
To validate the numerical solution, the DEP force was calculated and compared using both the Effective Dipole Moment (EDM) method and also by integrating the MST on the cell surface (radius = 5 µm) at different frequencies (Fig. 10). As shown, the results indicate the validation of the numerical solution and the inaccuracy of the EDM method. The AC electric field bears the advantage of a controlled DEP force due to its ability in tuning the applied frequency (Fig. 10).
The time average of electric stresses on the cell's surface is shown at the frequencies of 10 kHz and 10 MHz (Fig. 11). It can be concluded that the effect of a non-uniform AC electric field on the cell is different at various frequencies. Therefore, it is expected that the time average of the DEP force should also be different at these two frequencies.
The effects of buffer solution's dielectric properties on the magnitude of electric stresses at the frequencies 10 kHz and 10 MHz is shown in Figs. 12 and 13, respectively. From these figures, it can be seen that the magnitude of electric stresses increases as the ratio of ε ε / f p increases, and it decreases as the ratio of σ σ / f p increases. It can also be seen that the conductivity of the buffer solution has no effect on electric stresses at high frequencies.
The electric stresses from uniform electric field acting on the cell surface are shown in Fig. 14. In this figure, the electric field is considered to be in the y-direction. As can be seen, the resultant force exerted on a cell in a uniform electric field is zero, which is because the stress distribution over the cell surface is symmetrical.
Experimental results. Figure 15 shows the real part of the Clausius-Mossotti coefficient versus the applied electric field frequency to a red blood cell using the single-shell model with the parameters found in the literature 28,29 . The diagram was attained utilizing the spherical single-shell model defined by Gimsa et al. 30 .  www.nature.com/scientificreports www.nature.com/scientificreports/ In the previous section, the distribution of the electrical stresses on the cell surface was calculated. As shown (Fig. 14), these stresses are tensile and compressive at low and high frequencies, respectively. In this section, cell deformation in the presence of the electric field at various frequencies was investigated experimentally. Figure 16 shows cell deformation under the application of an AC electric field at different frequencies from 5 kHz to 20 MHz. According to Fig. 15, RBCs are in nDEP regime at the frequency of 50 kHz, while they are in pDEP regime at the frequencies of 500 kHz, 1 MHz, 3 MHz, 15 MHz and 20 MHz. As shown in Figs. 11 and 14, when a cell is affected in nDEP regime, the electrical stresses act on the cell surface in the form of tensile forces, whereas in pDEP regime, there will be compressive stresses exerted over the cell surface. According to Fig. 16, cell elongation increases as the frequency increases from 500 kHz to 1 MHz. The maximum cell deformation occurred at the frequency of 1 MHz, which corresponds to the peak of Clausius-Mossotti coefficient versus frequency plot (Fig. 15), while cell deformation was reduced at the frequencies of 15 MHz and 20 MHz. Also seen in Fig. 16, cells deformed under the effect of tensile and compressive electric stress distribution at the frequency of 50 kHz and 500 kHz, respectively. The cell deformation under tensile and compressive electrical stresses looks almost alike (in terms of elongation along y-axis). Looking at the tensile and compressive electrical stresses (Figs. 11 and 14), this result was expected. Based on these experimental results, no comprehensive discussion can be provided about the distinctive effects of tensile and compressive electrical distributions on the cell surface. Despite the electrical stresses magnitude, their tensile or compressive distribution on the cell surface may have an important effect on membrane structural changes. Further research is needed to address this issue.
Based on the above results, it can be concluded that the applied electric field frequency has a significant effect on cell deformation. To better show this effect, statistical analysis was conducted at different frequencies (from 300 kHz to 20 MHz) and the elongation ratio (i.e. the ratio of the cell elongation to its initial size) was calculated. The results from this analysis (Fig. 17), show that by increasing the applied electric field frequency up to 1 MHz, cell elongation ratio increases and reaches its maximum value. By further increasing the frequency, cell elongation ratio decreases. These observations match well with Clausius-Mossotti curve versus frequency plot (Fig. 15); by increasing the frequency from 300 kHz to 20 MHz the CM coefficient increases up to 1 MHz and reaches its maximum value and then decreases.

Conclusion
Cells are exposed to electromagnetic fields in various applications such as electroporation, electrofusion, etc. To better study these phenomena, it is necessary to understand the interaction between cells and electromagnetic fields. One of the most prominent electrokinetic effects of electric fields on cells is dielectrophoresis; a phenomenon that has found many biomedical applications. Previous works, mostly studied the dielectrophoretic effect of electric fields on a group of cells. Therefore, the present study focused on analyzing the behavior of single cells in an electric field. When a cell is exposed to an electric field, local changes in the electric potential distribution are induced. As a result, electrical stresses are exerted on the cell surface. These electrical stresses, in turn, lead to cell deformation and its membrane structural changes that trigger cell membrane permeabilization. Electrical stresses depend on the dielectric properties of cells, buffer solution, and the applied electric field characteristic. In the present study, the effect of applied voltage frequency on the distribution and the magnitude of electrical stresses exerted on a cell surface was numerically and experimentally studied. To study the interaction of cells and electric fields numerically, a novel numerical approach in solving the quasi-electrostatic equations in the presence of cells using the immersed interface method was applied. It was shown that the applied voltage frequency has a significant effect on cell deformation. Furthermore, it was shown that the effect of applied voltage frequency on cell deformation is consistent with the Clausius-Mossotti coefficient curve. In conclusion, a new approach to study the electroporation phenomenon was introduced. Developing an appropriate mechanical model for cell membranes can facilitate the studies conducted on the effects of electric field-induced cell deformations and membrane structural changes, which help to make biomedical processes and applications efficient. Moreover, the electrokinetic interaction of molecules and cells needs to be analyzed to study the transport of molecules across the cell membrane. The numerical approach applied in this research can help in analytical studies on cells and molecules interactions.
Ethics statement. In this research, healthy human red blood cells were used. Blood samples have been used in this research in an ethical manner. Also, authors confirm that all experiments were performed in accordance with relevant guidelines and regulations. Healthy human blood sample was obtained from a medical laboratory with the donor's written informed consent. The experimental procedure was approved by the research ethics committee, Pasteur Institute of Iran, reference number IR.PII.REC.1397.016.