Numerical study of electrohydrodynamic atomization considering liquid wetting and corona discharge effects Open Access

$A$

Area of the discharge surface, $m 2$

$C a e$

Dimensionless electric capillary number

$d j$

Jet diameter, $m$

$E$

Electric field, $−∇ϕ$

$E ∞$

Electric field intensity, $V/m$

$E o n$

Onset electric field

$E z$

Electric field component in the perpendicular direction to the outlet, $V/m$

$f e$

Electrostatic force, $N/ m 3$

$f σ$

Surface tension force, $N/ m 3$

$g$

Gravitational acceleration, $m/ s 2$

$I d$

Discharge current, $A$

$i.d$

Inner diameter, $mm$

$I cond$

Conductive current, $A$

$I conv$

Convective current, $A$

$i . d r$

Ring inner diameter, $mm$

$J$

Current density, $A/ m 2$

$n 0$

Uncorrected interface normal vector

$n c$

Corrected interface normal vector

$n w$

Boundary wall normal vector

$o.d$

Outer diameter, $mm$

$p$

Pressure, $Pa$

$Q$

Liquid low rate, $ml/h$

$r E$

Interelectrode distance, $mm$

$d fil$

Nozzle sharp edges fillet diameter, $mm$

$S$

Outlet area, $m 2$

$t$

Time, $s$

$u$

Fluid velocity, $m/s$

$u r$

Artificial compression term

$U z$

Liquid velocity component in the perpendicular direction to the outlet, $m/s$

$l c$

Cone length, mm

$v i$

Ionic wind velocity, m/s

$t j$

Time to first jet induction, s

$γ liq$

Phase fraction of liquid

$ε$

Permittivity, $F/m$

$ε 0$

Permittivity of free space

$ε r$

Dielectric constant

$θ$

Contact angle

$θ 0$

Angle between the uncorrected interface normal vector and the wall's normal vector

$κ$

Mean curvature of the free surface, $m − 1$

$κ e$

Electrical conductivity of the fluid, $S/m$

$µ$

Fluid viscosity, $Pa.s$

$μ e$

Mobility of charge of the gas

$ρ$

Fluid density, $kg/ m 3$

$ρ c$

Corona volumetric charge, $C/ m 3$

$ρ e$

Volumetric charge density, $C/ m 3$

$σ$

Surface tension, $N.m$

$ϕ$

Electric potential, $V$

$φ c$

Cone angle

$θ a$

Actual contact angle between the liquid's interface and the nozzle's wall

AC

Alternating current

CF

Carbon fiber

CSF

Continuum surface force

EHDA

Electrohydrodynamic atomization

N1

Nozzle 1

N2

Nozzle 2

VOF

Volume of fluid

e

Experiment

s

Simulation

Electrohydrodynamic atomization (EHDA), also named as electrospray, refers to a process where a fluid jet is broken up into particles/fibers (droplets) by a high electric field. Owing to its potential, electrospray has been applied in various areas, including pharmaceutics,^1–4 food industry,^5,6 fuel spray and combustion,^7,8 electric propulsion,^9,10 and jet printing.^11,12 Among the working regimes of electrospray ionization, the single cone-jet mode in which only one jet is induced is the most used due to its stability, controllability, and high yield rate. This mode was first investigated by Taylor¹³ who discussed the equilibrium condition of interfacial force with a specific semi-vertical angle (49.3°) of the conical shape, known as the Taylor cone. Notable early works include the publications of Rayleigh,¹⁴ Zeleny,¹⁵ and Taylor and Melcher^13,16 discussing fundamental aspects of electrospray. Taylor and Melcher¹⁶ suggested a mathematical model, named as the leaky-dielectric, to explain the relationship of physical parameters involved and approximate the electric force exerted on the liquid surface. This work has been utilized in various studies to describe the formation and development of the Taylor cone into droplets.

Due to limits of the analytical approach, numerical methods have been used to study the formation and behavior of electrospray.¹⁷ Several published works can be cited as follows: model the spraying of water and weak saline solutions at low-voltage using a double-layer nozzle;¹⁸ approximate the characteristics of a novel carbon fiber emitter electrospray for mass spectrometric study;¹⁹ and predict the electric potential in a large simulation domain.²⁰ Specifically, the numerical simulation of formation of the Taylor cone and its development into liquid jet and then droplets for various problems have attracted many researchers using different methods for different objectives. Herrada et al.²¹ proposed a new interface tracking method, and Du and Chaudhuri²² simulated the formation and breakup of droplets from axisymmetric charged liquid jets in electric fields. Several interesting physical phenomena occurring in charged droplet breakup and atomization, such as jet instability, necking, and the evolution of an unstable jet to droplet breakoff, were studied. Rahmanpour and Reza²³ and Panahi et al.²⁴ simulated the formation of the Taylor cone and into fluid jet under an external electric field²³ as well as the behavior of Newtonian viscoelastic jet.²⁴ Gamero-Castaño and Magnani²⁵ introduced numerical simulation employing Taylor's electric potential as a far-field boundary condition to achieve independence from geometry and potential of the electrodes, and the least-squared weighted residual method to the Young–Laplace equation for free surface tracking. Ouedraogo et al. found (i) the electrical properties of the liquid, such as electrical permittivity and conductivity can yield various mechanisms of fluid disintegration,²⁶ and (ii) the independence of the electrospray current on the applied voltage for conductive fluids.²⁷ Huh and Wirz²⁸ used a modified interpolation for electrical properties based on harmonic averages to improve the stability of the simulation of high conductivity fluids. Suo et al.²⁹ used a similar approach²⁸ to investigate the current transferred in the cone-jet mode for several different conductive liquids. These studies are most in good agreement with experimental results, strengthening the reliability of numerical simulations of the electrospray.

More interestingly, jet printing, one of the most common numerically investigated applications of electrospray, appears in numerous research works on various configurations, including conventional single nozzle^30–33 concentric non-conductive tip,³⁴ alternating current printer,³⁵ pulsating jet regime,^36,37 and coaxial nozzle.^38,39 All these works allow to study various physical behaviors of electrospray, which have enhanced the capability of electrospray in micro-/nano-manufacturing industry.

Previous studies concentrated on and validated the simulation of the Taylor cone by comparing its shape with experiment^{18,20,23,28,29,39} or focused on the development of the Taylor cone into jet under different conditions to explore the physics of electrospray.^{24,31,36–38} Experimental observations show that during electrospray, liquid tends to expand and cover the capillary nozzle's tip due to liquid wetting. However, so far numerical studies on the electrospray have either briefly mentioned or ignored the phenomenon, yielding inexact simulation of the cone shape and dimensions and inaccurate estimation of surface tension force based on the liquid's interface curvature and nozzle diameter. Furthermore, since corona discharge from nozzle electrode significantly affects the formation of the Taylor cone and the steady cone-jet mode, the introduction of corona discharge into simulations allows us to possibly approach real-life electrospray. Therefore, in this work, we perform electrospray simulations with a new approach in which the liquid wetting and corona discharge are additionally examined. The liquid wetting, considered as one of the key parameters, allows us to improve the simulating accuracy of the Taylor cone shape as well as its formulation and development. Meanwhile, electrospray–corona discharge integration is conducted to investigate the effects of corona discharge on electrospray and consolidate the importance of liquid wetting in electrospray analyses. The accuracy and reliability of the present simulation are investigated and evaluated by published works and, especially, experimental works using a high-speed camera to capture the formation and development of the Taylor cone.

The numerical simulation for the jet formation at the nozzle tip is investigated and evaluated by experiment using the nozzle-to-ring configuration whose schema is described in Fig. 1. Liquid is fed to the tip of a nozzle, where its interface is deformed by an applied electric field before elongating into a jet. The system consists of a steel nozzle, a metal ring with an outer diameter of 10 mm and an inner diameter of 4.8 mm, a high-voltage power source, and a syringe pump. Two types of nozzles, Nozzle 1 (N1) and Nozzle 2 (N2), are employed for the experiment, with N1 featuring outer diameter (o.d) and inner diameter (i.d) of 0.3 and 0.15 mm, respectively, and N2 featuring o.d and i.d of 0.7 and 0.4 mm, respectively (stainless-steel G29 nozzle and G22 nozzle—Musashi Engineering, Inc., Japan). The distance between nozzle and ring (⁠ $r E$⁠) is adjusted at 4 mm for N1 and 7 mm for N2 to establish sufficient electric field. A high-speed camera (HHC X9 Pro, Mega Speed, USA) and a microscope camera (AM 4113, Dino-Lite, Australia) are used to capture the formation of the Taylor cone and its development into fluid jet and droplets at the nozzle tip. The used liquid is polyethylene glycol with a molecular weight of 200 (PEG-200 - Merck, Australia).

The objective of this experiment is to observe the formation of Taylor cones, jet elongation, and propagation of liquid using both high-speed and microscope cameras. To achieve a stable Taylor cone and jet, the applied voltage is fixed at 3400 and 5400 V for N1 and N2, respectively. Meanwhile, a voltage of 5600 V is applied to capture the formation process of a Taylor cone with the high-speed camera. The flow rate is fixed at $Q=0.75 ml/h$⁠, which is determined by the approximation of the minimum flow rate required for a stable Taylor cone $Q∼σ ε r ε 0/ρ κ e.$⁴⁰ These experimental conditions allow for the observation and analysis of the physical phenomena of interest. Moreover, the contact angle of PEG-200 is determined by capturing the droplet on a stainless-steel plate with a digital camera (HHC X9 Pro, Mega Speed, USA). The measurements are performed using the ImageJ software (NIH). Here, similar material with metal nozzles is employed to reflect equivalent conditions with the liquid on nozzles.

We utilize the interFOAM solver of the open-source OpenFOAM package. This solver has been proven reliable in a wide range of complicated problems, such as cavitation or microfluidic capillary flows,⁵⁵ and, hence, emerges as the suitable solution for electrospray simulation. The additionally developed electrostatic solver and the existing fluid dynamic solver are coupled and solved consecutively each time step. A detailed description of the solving procedure is presented in Fig. 2.

Simulations in this study involve two nozzle configurations, Nozzle 1 (N1) and Nozzle 2 (N2), whose sizes and conditions conform to experimental setups. Figure 3 shows the computational mesh and dimensional annotations of these configurations. The feeding nozzle (i) is positive electrode with constant voltage $ϕ= ϕ 0$⁠. The ring electrode (ii − i.d_r = 4.8 mm) is set with $ϕ=0 V$⁠. The domain is enclosed by atmospheric boundaries with far-field conditions (more information of the computational model is in supplementary material S4.1 and S4.2). An axisymmetric 2D model is employed to simplify the solution and reduce computational resources, with the front and back patches having a wedge condition with a 5° angle.

The computational mesh and domain, shown in Fig. 3, is constructed by the open-source gMsh mesh generator developed by Geuzaine and Remacle.⁵⁶ In this work, hybrid mesh resolution is employed to reduce the number of mesh elements whilst maintaining good flow alignment. The nozzle's sharp edges are filleted to guarantee better shape agreement with empirical nozzles as well as to describe somewhat more accurately the motion of the fluid on the surface of the capillary nozzle (see supplementary material S2.2). Here, total numbers of cells of 49 620 and 68 551 for N1 and N2 are obtained as minimum cell size is kept at $3 μm$ in both cases. For numerical schemes, we employ Gauss minmod for divergence components, Gauss linear corrected for Laplacian components, and cellLimited Gauss linear for gradient components. These schemes seem to provide good stability with an average time step of $3× 10 − 7 s$⁠. Moreover, the physical properties of air and PEG-200 used in this work are shown in Table I.

The contact angle for nozzle is measured by experiment and used as the boundary condition for contact angle in simulations. Figure 4(a) shows (i) the experimental capture of a stationary PEG-200 droplet on a stainless-steel plate and (ii)(vi) verification simulation results of different contact angle conditions. Comparing each simulation case to experimental manifestation specifies most comparable shape between experiment and (iii) $θ s=20°$⁠. On top of that, measurement from the experimental capture also yields $θ e=20°$⁠, potentially indicating the reliability of the contact angle correction discussed in the previous section. Furthermore, the influence of the contact angle on the propagation mechanism of liquids is evaluated by fluid-dynamic (no electrostatic field) simulation results provided in Fig. 4(b). These results are obtained from completely similar boundary conditions, except for the contact angle, and at the same time instant. Case with no contact angle considered yields a liquid's surface being attached to the inner edge of the nozzle's tip, as commonly seen in preceding literature. Meanwhile, different contact angles result in a different wet tip's surface. Particularly, the smaller the liquid's contact angle is, the quicker its surface would propagate on the nozzle's wall and reach its outer edge. Therefore, to ensure adequate agreement with empirical state, the condition $θ=20°$ is applied on the nozzle boundary in subsequent electrospray simulations.

Corona discharge is considered as the inherent process of electrospray and may take various forms.⁵⁹ Discharge can initially originate from the metal electrode,⁶⁰ the tip of the cone where electric field intensity is high,⁶¹ or molecular ions emission from jet/droplets breakup,^62,63 and can be measured in both dry-spray (no liquid flow rate) and wet-spray conditions.^59,64 Since the formulation process of the Taylor cone from a dry nozzle (no liquid on the nozzle's tip) is our concern here, we would assume that the corona discharge in our case solely originates from the outer sharp edge of the metal nozzle. The outer sharp edge is the closest point of the nozzle to the ring electrode, causing the electric field to be most intense at this position. Furthermore, the sharpness of the edge, despite filleted, can emphasize electric field intensity in the adjacent region. Given that corona discharge is a consequence of strong electric field and the feasibility of discharging from metal surface, we consider corona discharge from the outer edge of the nozzle in our analysis.

The empirical and numerical cone shapes for the two nozzle configurations are digitalized and plotted in Fig. 5. It is worth noting that in our simulation, the voltage in case N2 has to be $ϕ s=6500 V$⁠, approximately 20% higher than in experiment (⁠ $ϕ e=5400 V)$⁠, to maintain steady cone-jet mode. It could be the 2D axisymmetric simplification of our model and the errors in experimental electrode arrangement that cause this difference. Nonetheless, our results show that the model that additionally considers liquid wetting and contact angle has produced visually comparable shapes with experiment in both cases. For the N1 configuration, experiment cone length $( l c)$ is at around 0.162 mm, while simulation registers a cone length of 0.142 mm (12% difference). This value for the N2 configuration is 0.45 and 0.41 mm for experiment and simulation, respectively, marking a 9% deviation. Simultaneously, the cone angles (⁠ $φ c$⁠) of the two sets of data also show good agreement in both cases. For N1, cone angles are 82° in experiment and 91° in simulation, while in case N2, we obtain 84° and 94°, respectively. Moreover, empirical images showed that the jet produced in the N1 configuration is visibly thicker than in the N2 configuration, which can also be observed in simulation, indicating good qualitative agreement in terms of jet diameter between the results. Quantitatively, our measurement from the digitalized cones yields experimental jet diameter $d je=13.8 μm$ and numerical jet diameter $d js=7.96 μm$ for the N1 configuration, while these parameters are $d je=14.2$ and $d js=9.8 μm$ for the N2 configuration. Overall, the inclusion of the contact angle and subsequent liquid propagation on the nozzle's wall accurately recreated the liquid state empirically manifested despite some deviations from experiment. Errors in voltage and measured quantities or inability of the model to produce the jet disintegration into a spray cloud as in experiments may be due to the axisymmetric model, grid resolution, and inaccuracy in image processing. These defects are to be addressed in future works. Considering the proposed scope of this study and the presented correlations between our experiments and simulations, the developed computational model can be deemed satisfactory and employed in next investigations.

The correlations between current and jet diameter with conductivity, voltage, surface tension, flow rate, and viscosity are widely discussed operational characteristics of electrospray. To optimize computing resources, we limit our simulations to employing solely the Nozzle 1 configuration, which enables us to investigate parametric relationships and confirm the dependability of our recently developed solver.

The correlation between spray current and jet diameter with electrical conductivity $κ e$ is graphed in Fig. 6. Simulations yields $I∝ κ e 0.22$ and $d j∝ κ e − 0.18$ in the range $0.01 mS/m≤ κ e≤0.5 mS/m$⁠, which agrees well with early experimental reports by Tang and Kebarle⁶⁵ (⁠ $I∝ κ e 0.2 − 0.4$⁠) for $κ e> 10 − 2 mS/m ( 10 − 7 Ω − 1 cm − 1)$⁠, and does not divert much from the asymptotic law by Gañán-Calvo⁶⁶ (⁠ $I∝ κ e 0.5$⁠; $d j∝ κ e − 0.17$⁠). Qualitatively, this can be a consequence of the thicker charge concentration because of high conductivities, which exerts more Coulombic force under the influence of electric field and, subsequently, gives rise to higher jet velocities. The higher flow rate, coupled with the conservation of mass, leads to a reduction in the diameter of the jet, resulting in the commonly observed electrical conductivity effects. (Additional simulation results are provided in supplementary material S4.3.)

Furthermore, the impacts of applied voltage and the liquid's surface tension are also comparable with other reports. Our model predicts an independence of current on the applied voltage and a reduction of jet diameter in the steady cone-jet range, as shown in Fig. 7(a). Here, simulations yield a constantly developing Taylor cone shape for lower voltage range (⁠ $ϕ<3400 V$⁠) due to lack of equilibrium between feeding and spraying flow rate. We consider this range as unsteady cone-jet and neglect it in any dimensionless analyses. Under steady conditions, simulations have successfully reproduced the independent current–voltage relationship as has been reported in recently published experiment data.^67,68 The slight decrease in jet diameter with increasing applied voltage also shows good conformity with experimental gains of Gañán-Calvo et al.⁶⁹ While the jet diameter characteristic is caused by a growth in electrostatic force, the voltage–current relation results from the compensation of the rise jet velocity for the contraction in diameter. As suggested by Eq. (19), an increase in velocity and a corresponding decrease in diameter would result in an unchanged current (see additional results in supplementary material S4.3). It is worth mentioning that our model also produces charge density $ρ e$ reliance with the phase fraction value $γ liq$⁠, which is used to calculate jet diameter; hence, both the conductive and convective currents are affected by jet diameter reduction. For the surface tension effects in the cone-jet mode, scaling law from⁶⁶ suggests analogous trends with conductivity (⁠ $I∝ σ 0.5$⁠; $d j∝ σ − 0.17$⁠). As can be seen in Fig. 7(b), in the steadily working range, our model registers mildly consistent $I∝ σ 0.34$ for current, while jet diameter shows an agreed downward trend.

The spray current and jet diameter variation with flow rate and viscosity are shown is Fig. 9. Our simulations predict a similar upward trend $I∝ Q 0.92$ and $d j∝ Q 0.2$ in the steady cone-jet range, which agrees in order with scaling law⁶⁶ (⁠ $I∝ Q 0.5$⁠; $d j∝ Q 0.5$⁠), yet overpredict the change in current and underestimate the increase in jet diameter. On the one hand, numerical results show negligible change in jet diameter as viscosity increases, which is comparable with reported work⁷² predicting the insignificant effect of viscosity. Furthermore, scaling law from Gañán-Calvo⁶⁶ excludes viscosity on both spray current and jet diameter relations, implying the independence of these parameters on liquid's viscosity. On the other hand, our spray current results register an upward trend as viscosity value rises, which might have been mildly overestimated considering its insensitivity to viscosity.^66,73 This dissimilarity between simulations and empirical data on spray current was also obtained by Suo et al.²⁹ in which surface charge convection and strong tangential stress along the cone-jet surface were mentioned as the cause of the increasing trend of current. This explanation can also be applied in our study where we record a reduction of surficial charge density as viscosity increases due to the weakened recirculating motion (see supplementary material S4.3). This motion arises due to the tangential electric stress along the Taylor cone's surface and is weakened or canceled when the liquid's viscosity is higher, which is consistent with experimental works of Barrero et al.^74,75

To achieve better insight into the influence of contact angle on the cone-jet mode's characteristics, simulations on three contact angles conditions $θ=10°,20°, and 40°$ have been performed for both the N1 and N2 nozzle with $ϕ=3400 V$ and $Q=0.75 ml/h$⁠. Figures 10(a) and 10(b) present Taylor cone's shape of N1 and N2 configuration, respectively, whereas Fig. 10(c) provides spray current and jet diameter variation with contact angles. Our results obtain negligible differences in the Taylor cone's shape despite considerable changes in liquid's contact angle, which results in the corresponding somewhat unaffected spray current and jet diameter characteristics. In both nozzles, our simulations produce relatively similar cone's shapes for all three contact angle conditions, except for slight differences at the attachment positions at the filleted outer edge of the N1 setting [Fig. 10(a)]. However, this effect only localizes to liquid region adjacent to the nozzle's wall and does not affect the rest of the Taylor cone. Moreover, the actual angles between the liquid's surface and the wall of the nozzle's tip $θ a$ [annotated in Fig. 10(b)] are considerably higher than the contact angle conditions employed and even approaches $90°$ in some cases. Underlying reasons for these behaviors can be the overwhelming effects of other factors, namely applied voltage, flow rate, and other liquid's properties (conductivity, viscosity, and density)⁷⁶ on the cone's shape. For a chosen liquid, this potentially indicates that the roles of contact angle can be promoted under different voltage and flow rate. To avoid further complicating the investigation, we only conduct contact angle analyses under two separated conditions, one with a higher voltage (⁠ $ϕ=3600 V$⁠; $Q=0.75 ml/h$⁠) and the other a lower flow rate (⁠ $Q=0.6 ml/h;ϕ=3400 V$⁠), using the N1 nozzle.

As observed in experiments,^77–80 higher applied voltage and lower flow rate can shorten the length of the Taylor cone $l c$⁠, enlarge the cone angle $φ c$⁠, and reduce the actual contact angle $θ a$⁠. Such behaviors of cone length and cone angle occur to balance the increasing electric stress due to rising voltage as well as the decreasing feeding volume due to lessening flow rate, which eventually keeps the cone-jet mode stabilized.⁷⁹ 

Figure 10(d) shows spray current and jet diameter variation with contact angle for new voltage and flow rate condition while the cone's shape for each separated condition is provided in Figs. 10(e) and 10(f), respectively. One the one hand, in both voltage and flow rate conditions, no significant difference in spray current, and jet diameter as well as cone's shape is obtained between $θ=10°$ and $20°$⁠. At these contact angles, our simulations produce noticeably shorter $l c$⁠, larger $φ c$ compared to cases in Figs. 10(a) and 10(b), conforming with the previously discussed experimental observations. Here, the actual contact angle $θ a$ is also smaller due to higher voltage and lower flow rate, yet still significantly larger than the declared contact angle of liquid (⁠ $θ=10°$⁠, $20°)$⁠. On the other hand, in case $θ=40°$⁠, the higher voltage and lower flow rate both yield receding cones with markedly higher spray currents and larger jet diameter [Fig. 10(d)]. At $θ=40°$⁠, the surface tension induces force more inclined toward the center of the cone, whereas at $θ=10°$ and $20°$⁠, this force is diverted more toward the nozzle's wall. Therefore, instead of producing a shortened cone, the high contact angle $θ=40°$ develops a receding motion of the liquid's interface on the nozzle in which the liquid's surface quickly withdraws from the attachment position at the filleted outer edge, shrinking the Taylor cone radially. Despite not commonly reported, this retractive response can still be encountered in experiments.⁸¹ 

Overall, our results could imply that although the influence of other factors, such as voltage and flow rate, is more decisive in terms of modifying the Taylor cone's shape and the characteristics of the cone-jet mode, the liquid's contact angle is still impactful in certain circumstances, i.e., high applied voltages and low flow rates. The quantitative correlation between the receding mechanism with the liquid's contact angle $θ$ as well as the actual contact angle $θ a$ remains an open problem and can be explored in dedicated research.

Figure 11 shows different liquid states at six consecutive time steps obtained from our experiments and simulations for the N2 configuration. Figure 11(a) includes numerical results neglecting any corona calculation, while Fig. 11(b) provides results considering corona discharge with a current of $I d=1 μA$⁠. By further considering corona discharge, our simulations show better agreement with empirical captures displayed in Fig. 11(c). As stated, experimental data yield a slowly propagating liquid [at t₁ and t₂ in Fig. 11(c)] on the nozzle and the deformation of liquid's interface into steady cone-jet mode can only be obtained after the liquid completely covers the outer edge of the nozzle (t₃ and t₄). Similar behavior is reproduced by electrospray–corona simulation [Fig. 11(b)] where a cluster of analogously signed charge can be observed forming around the electrode (see Fig. 12), reducing the field intensity on the liquid surface to below its critical value for jet induction. As the nozzle edge is progressively covered by propagating liquid, the surface discharge is, therefore, obstructed and alleviated, causing an increase in electric field intensity. Due to this relationship, critical field value can be reached at t₄, gradually deforming the interface, and giving rise to jet afterward [at t₅ and t₆ in Figs. 11(b) and 11(c)]. On the other hand, if corona discharge is ignored, the field strength induces jet at very early time steps [t₂ to t₄, Fig. 11(a)] when the liquid is still expanding on the nozzle, which is not consistent with empirical records.

In our simulations, we have shown that corona discharge has a delay effect to the formation of the cone-jet mode by weakening the electric field. Based on this, we reasonably predict that if the discharge is vigorous enough or originates from other positions on the nozzle (due to microscale surficial imperfections), steady cone-jet mode may not be achieved. Moreover, in cases where corona discharge and the cone-jet mode can coexist, further increase in applied voltage would not destabilize the cone-jet mode also because of the electric field reduction effect. These predictions are on par with many electrospray–corona research works,⁵⁹ indicating good physical conformity with empirical gains of our numerical model.

Figure 12 provides additional results on the velocity field and the corona discharge contour at two greatly distinguished time instants. The discharge gives rise to a large cluster of charge (charge cloud) and an ionic wind field at t₂. Ionic wind is the movement of air caused by the momentum transferred from electrically charged molecules to neutral ones when they are accelerated by an electric field,⁸² and it is the commonly known by-product of corona discharge. In this work, the average velocity of the ionic wind measured at the outlet boundary stays at around $v i=2.6 m/s$ for unobstructed discharge (t₂ in Fig. 12), which concurs in order with empirical measurements of pure corona-induced ionic wind of Zhang et al. (needle-to-ring: 1–8 m/s⁸³) and Wang et al. (multiple-needles-to-mesh: 0.7–1.7 m/s⁸⁴ and two non-uniform needles-to-ring: 0.75–3.5 m/s⁸⁵) or of our team's recently published works (dual pin discharge: 0–1.3 m/s).^47,49 At t₅, the discharge coverage and the velocity are significantly lessened, influencing the system in the manner described earlier. Additionally, the inset in Fig. 12 shows a long exposure capture of the corona discharge in the case of a dry nozzle, showing a glowing ionized region nearby and underneath the outer edge. This glow is the light emitted by excited neutral atoms when their electrons revert to ground state and can behave differently under different surrounding gaseous fields.⁸⁶ The characteristic glow indicates strong electric field, verifying the intense field in the vicinity of the relatively rounded edge in the contour results of Figs. 11(a) and 11(b) and confirming our earlier assumption regarding the origin of corona discharge in the electrospray.

The electrospray–corona simulation model developed in this study, with the proposed assumptions, is able to accurately reproduce the empirical formation of the Taylor cone and achieve reasonable consistency with existing literature regarding the effects and behaviors of corona discharge and ionic wind. The inclusion of corona discharge in the electrospray simulation highlights the crucial role of liquid wettability, which should be taken as a prominent impacting factor. Neglecting this factor may cause the solution flawed and unable to achieve a steady cone-jet mode.

Additionally, we perform multiple simulations with different discharge currents to further explore the effects of corona to the formulation of the Taylor cone. Figure 13(a) plots the ionic wind velocity and the maximum electric field (max E-field) at the nozzle's tip at $t m=30 ms$ when the liquid has not yet reached the outer edge of the nozzle. Results for discharge-free case (⁠ $I d=0 μA$⁠) are also included to provide a thorough comparative view. Simulations yield that both ionic wind velocity and max E-field follow a similar exponential relationship with discharge current, yet opposite in trend (⁠ $v i∝ I d 0.34$⁠; $E max∝ I d − 0.54$⁠). Simulated wind velocity is within 1–3 m/s, which still registers good order agreement with the discussed literature. In terms of electric field, simulations predict the field within 0.5 × 10⁷ – 2 × 10⁷ V/m, which is weakened from $E max$ ≈ 2.5 × 10⁷ V/m in the discharge-free case. Here, we also calculate the time required for first jet induction t_j (elapsed time from the initial instance to the beginning of the induction of the first jet) to investigate the influence of discharge current on the delay effect of corona. At low current $I d=0.1 μA$⁠, the weakened field can still give rise to jet when the liquid is still propagating on the nozzle's surface, while in the range $0.5 μ A ≤ I d≤1.5 μA$⁠, the propagation of liquid and the delay effect can be clearly observed. In this range, the time t_j is relatively close with only several milliseconds difference (52; 57; 65; and 71 ms). This tendency suggests that higher discharge current requires a larger portion of the discharge area to be covered by liquid for electric field to sufficiently increase and exceed the critical value for jet induction.

Moreover, Fig. 13(b) provides a graph of cloud charge for three distinguished discharge currents, (i) $I d=0.1$⁠, (ii) $I d=0.5$⁠, and (iii) $I d=1.5 μA$ also at the same time step $t=30 ms$⁠. Broader cloud width is generated in cases of higher corona discharge current. As the discharge current increases, the Coulombic force in the charge cloud is more intense, expelling the charges outward, resulting in wider charge clouds (additional contour results in supplementary material S4.3).

In our study, the current values at which delay effect is observed (⁠ $0.5 μ A ≤ I d≤1.5 μA$⁠) are weak in intensity, yet realistically possible considering corona discharge current may vary from a few hundred nanoampere to a few dozen microampere for various configurations.^{45,46,48,49,51,59} Moreover, the 2D axisymmetric model employed in this work may overemphasize the field reduction effect of the charge cloud from corona and the ionic wind velocity; thus, to avoid greater deviation from reality because of overprediction of these factors, low current range is investigated. Appropriate 3D models can be utilized in the future for more accurate numerical reproductions of corona discharge in electrospray.

A relevant numerical simulation for the cone-jet mode of electrospray has been carried out using a new approach in which the liquid wetting and corona discharge are considered as impacting parameters. Numerical results show that the present simulation achieves Taylor cone shapes very similar to those by experiments of two different configurations. Also, the electrospray's characteristics by the new approach, including the relations of the spray current and jet diameter with respect to the electrical conductivity, applied voltage, surface tension, flow rate, and viscosity are in reasonable agreement with those stated in published works.^{65–69,72,73} The impacts of contact angles on the cone-jet mode are shown to be minor compared to those of applied voltage, flow rate, and other properties of liquids and, however, can be significant under higher voltages and lower flow rates. Furthermore, the electrospray–corona discharge integrated simulation, evidenced by experimental results, depicts that the introduction of corona discharge yields better and more appropriate simulation of the Taylor cone formation with higher accuracy compared with one by conventional non-corona simulations. The work demonstrates the crucial role of liquid wetting and corona discharge not only in the cone shape formation of the electrospray but also in further improving the numerical study of electrospray in the relation with corona discharge.

See the supplementary material for more information about the experiment installation, the contact angle correction formulation, additional experiment results, the computational model, jet diameter estimation, and more simulation results.

We acknowledge the support of time and facilities from Ho Chi Minh City University of Technology (HCMUT), VNU-HCM for this study. This work is a collaboration between VNU-HCM Key Lab. for Internal Combustion Engine, Ho Chi Minh City University of Technology, Ho Chi Minh City, Vietnam and School of Engineering and Built Environment, Griffith University, Queensland, Australia.

The authors have no conflicts to disclose.

The data that support the findings of this study are available from the corresponding authors upon reasonable request.