2 Experimental characterization

Figure 2.

(a) Sketch of the experimental set-up. (b) Typical current versus flow rate measurements, down to the minimum stable flow rate, for a propylene carbonate solution with $K=0.0282~\text{S}~\text{m}^{-1}$ ; each series is taken at constant emitter potential. (c,d) Geometry of the cone jet near the minimum flow rate for emitter potentials of 3400 V and 4000 V.

Figure 2 illustrates the experimental set-up and the methodology for measuring minimum flow rates. The liquid is stored in a reservoir, and fed to the emission site through two fused silica tubes connected in series by means of a zero-volume interface. The free end of the second tube is made electrically conducting with sputtered iridium and, upon applying a potential difference $V_{E}$ with respect to a facing collector disk, anchors the base of the cone jet, see figure 2(c,d). The outer diameter of this tube, $360~\unicode[STIX]{x03BC}\text{m}$ in all experiments, is also the diameter of the base of the cone jet. The flow rate is controlled by applying a known pressure to the reservoir, while the emission site is open to the atmosphere. For accurate flow rate measurement this pressure difference $\unicode[STIX]{x0394}P$ needs to be sufficiently high throughout the range of flow rates, which covers many orders of magnitude for the solutions investigated. To accommodate this we take advantage of the strong dependence of the flow rate on the inner diameter of the tube, $Q\sim D^{-4}$ , and select capillary tubes with different inner diameters depending on the solution that is electrosprayed. The fused silica tube with its free end inserted in the fluid reservoir has an inner diameter of $230~\unicode[STIX]{x03BC}\text{m}$ , its contribution to the overall hydraulic resistance is insignificant, and works as a flowmeter: to obtain the $Q(\unicode[STIX]{x0394}P)$ calibration of the line we insert air and measure the velocity at which the liquid–air interface moves at fixed $\unicode[STIX]{x0394}P$ . The electrospray current is measured at the collector with an electrometer, and the pressure of the reservoir with an MKS Baratron pressure transducer. Their analogue outputs are recorded with a data acquisition card, a computer and a LabView application. Figure 2(b) illustrates the characterization of the minimum flow rate for a solution of propylene carbonate with a conductivity of $0.0282~\text{S}~\text{m}^{-1}$ : in a typical run the emitter voltage is fixed and the reservoir pressurized to make the flow rate a few times larger than the minimum; the pressure is then slowly released through a needle valve, while recording it together with the current (see inset in figure 2 b); the recording is stopped once the electrosprays becomes unstable, as identified by the collapse of the electrospray current; these steps are repeated at several emitter potentials, and at least twice for each potential. Figure 2(b) shows current versus flow rate curves near the minimum flow rate for eight independent runs, taken at four different emitter potentials. Square symbols are indicative of the point at which the current collapses. The average of the eight minimum flow rates is $3.10\times 10^{-12}~\text{m}^{3}~\text{s}^{-1}$ and the relative standard deviation is 4.3 %. This narrow dispersion suggests the average value to be a deterministic output. In particular, there is no trend between the minimum flow rate of each run and the potential of the emitter. Note also that, as illustrated by figures 2(c) and 2(d), the geometry and volume of the meniscus change significantly with the emitter potential: near the instability condition the meniscus at $V_{E}=3400$  V remains symmetric, conical, and with an aspect ratio near one, while the meniscus at $V_{E}=4000$  V is highly stressed, flatter, and its apex is far from the axis of the emitter. The independency between the minimum flow rate and the geometry and size of the meniscus is to be expected from the large ratio between the radii of the emitter and the jet. The characteristic radius of the jet at the minimum flow rate is $r_{c}=0.080~\unicode[STIX]{x03BC}\text{m}$ , i.e. the ratio between the radii is 2250. Under these conditions the minimum flow rate is a function of the physical properties of the fluid only, and does not depend on the volume of the meniscus or the potential of the emitter.

We have measured the minimum flow rate of solutions of tributyl phosphate (TBP) and propylene carbonate (PC). The values of the density, surface tension, viscosity and dielectric constant at $25\,^{\circ }\text{C}$ are $976~\text{kg}~\text{m}^{-3}$ , $0.0271~\text{N}~\text{m}^{-1}$ , 0.00339 Pa s and 8.91 for tributyl phosphate, and $1195~\text{kg}~\text{m}^{-3}$ , $0.0414~\text{N}~\text{m}^{-1}$ , 0.00253 Pa s and 64.9 for propylene carbonate (Riddick, Bunger & Sakano Reference Riddick, Bunger and Sakano1986). The very different dielectric constants make it possible to evaluate the effect of this parameter on the minimum flow rate. Table 1 lists the electrical conductivity of each solution (the pure liquids are doped with an ionic liquid to obtain the desired conductivity), together with the Reynolds number, the dimensionless minimum flow rate $\unicode[STIX]{x1D6F1}_{min}$ and the characteristic radius of the base of the jet $r_{c,min}$ at this point. The first obvious observation is that the minimum dimensionless flow rate for both liquids increases significantly with the electrical conductivity, contrary to the predictions of every criterion except for (1.4). Figure 3 plots the minimum flow rate versus the Reynolds number. The height of the error bar is twice the standard deviation of the measurements. The minimum flow rate increases at decreasing Reynolds number for both liquids. At fixed Reynolds number propylene carbonate exhibits a significantly larger minimum flow rate than tributyl phosphate. The experimental data are incompatible with laws (1.1) and (1.2), which do not depend on the Reynolds number. The experimental data used by Fernández de la Mora & Loscertales (Reference Fernández de la Mora and Loscertales1994) and Gañán-Calvo et al. (Reference Gañán-Calvo, Dávila and Barrero1997) to justify (1.1) and (1.2) only support them as approximate criteria. For example figure 9 in Fernández de la Mora & Loscertales (Reference Fernández de la Mora and Loscertales1994) shows measurements for water at varying conductivity suggesting that $\unicode[STIX]{x1D6F1}_{min}$ increases at decreasing Reynolds number, similarly to what we observe in our experiments. Figure 9 in Chen, Pui & Kaufman (Reference Chen, Pui and Kaufman1995), also for water solutions at varying conductivity, suggests that $\unicode[STIX]{x1D6F1}_{min}$ is constant for water. However the Reynolds numbers for these solutions are between 1.28 and 3.4, i.e. they fall in the flatter section of the $\unicode[STIX]{x1D6F1}_{min}(Re)$ curves shown in figure 3 where $\unicode[STIX]{x1D6F1}_{min}$ does not change substantially. Note also that $\unicode[STIX]{x1D700}Re$ in figure 3 ranges between 16.0 and 0.654 for tributyl phosphate, and between 200.9 and 22.1 for propylene carbonate. Thus, the minimum flow rate criteria proposed by Gañán-Calvo et al. (Reference Gañán-Calvo, Rebollo-Muñóz and Montanero2013) and Higuera (Reference Higuera2017) seem unable to explain the experimental relationship because most investigated solutions have large $\unicode[STIX]{x1D700}Re$ values and would need to follow the Reynolds-independent laws (1.3) and (1.5) respectively. The disagreement is especially prominent for the propylene carbonate solutions due to their very large $\unicode[STIX]{x1D700}Re$ : for these very large values of $\unicode[STIX]{x1D700}Re$ the viscous stress should be much smaller than the polarization force in the case of Gañán-Calvo’s criterion, much smaller than inertia in Higuera’s criterion, and the Reynolds number should drop from the $\unicode[STIX]{x1D6F1}_{min}(Re,\unicode[STIX]{x1D700})$ relationship contrary to the experimental data in figure 3. A possible explanation of the apparent disagreement with Gañán-Calvo et al. (Reference Gañán-Calvo, Rebollo-Muñóz and Montanero2013) and Higuera (Reference Higuera2017) is that the $\unicode[STIX]{x1D700}Re$ values investigated may not be large enough for the inviscid limits (1.3) and (1.5) to apply. Furthermore, the Reynolds number for these solutions (and probably for most electrosprayable solutions) is never too large or too small, which limits the accuracy that can be expected from asymptotic estimations for very large or very small Reynolds numbers. Finally, we note that a second important observation from figure 3 is that the dependency on the Reynolds number requires the viscous stress to play an important role in the determination of the minimum flow rate.

Figure 3.

Dimensionless flow rate versus Reynolds number for solutions of propylene carbonate and tributyl phosphate. The height of the error bar is twice the standard deviation of the measurement.

Table 1.

Electrical conductivity, Reynolds number, dimensionless minimum flow rate and characteristic radius of the base of the jet of the propylene carbonate (PC) and tributyl phosphate (TBP) solutions characterized in experiments.

3 Numerical model

We next summarize the leaky dielectric model particularized to cone jets. A detailed description of the model and the numerical algorithm can be found in Gamero-Castaño & Magnani (Reference Gamero-Castaño and Magnani2019). The variables of the model are the velocity $\boldsymbol{v}(\boldsymbol{x})$ and pressure $p(\boldsymbol{x})$ fields in the liquid, the position of the free surface $R(x)$ , the surface charge $\unicode[STIX]{x1D70E}(x)$ and the potential fields inside $\unicode[STIX]{x1D719}^{i}(\boldsymbol{x})$ and outside $\unicode[STIX]{x1D719}^{o}(\boldsymbol{x})$ the liquid. The electric field is given by $\boldsymbol{E}=-\unicode[STIX]{x1D735}\unicode[STIX]{x1D719}$ . The equations are expressed in dimensionless form using the scales $r_{c}$ , $v_{c}=Q/(\unicode[STIX]{x03C0}r_{c}^{2})$ , $p_{c}=\unicode[STIX]{x1D6FE}/r_{c}$ , $I_{c}=(\unicode[STIX]{x1D6FE}KQ)^{1/2}$ , $E_{c}=I_{c}/(\unicode[STIX]{x03C0}r_{c}^{2}K)$ , $\unicode[STIX]{x1D719}_{c}=r_{c}E_{c}$ and $\unicode[STIX]{x1D70E}_{c}=I_{c}/(2\unicode[STIX]{x03C0}r_{c}v_{c})$ . Henceforth symbols capped with a tilde will designate dimensional variables. The leaky dielectric model with constant fluid properties includes the equations of conservation of mass, Navier–Stokes and conservation of charge at the surface and the Laplace equation for the electrical potential (Melcher & Taylor Reference Melcher and Taylor1969; Saville Reference Saville1997; Higuera Reference Higuera2003):

(3.1) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{v}=0,\quad \boldsymbol{v}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{v}=-\frac{\unicode[STIX]{x03C0}^{2}}{\unicode[STIX]{x1D6F1}^{1/2}}\unicode[STIX]{x1D735}p+\frac{\unicode[STIX]{x03C0}}{\unicode[STIX]{x1D6F1}^{1/2}Re}\unicode[STIX]{x1D6FB}^{2}\boldsymbol{v}, & \displaystyle\end{eqnarray}$$

(3.2) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{d}}{\text{d}x}(Rv_{s}\unicode[STIX]{x1D70E})=2RE_{n}^{i}(1+R^{\prime 2})^{1/2}, & \displaystyle\end{eqnarray}$$

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D719}^{i,o}=0, & \displaystyle\end{eqnarray}$$

where $v_{s}$ stands for the velocity at the surface, $E_{n}$ and $E_{t}$ for the components of the electric field normal and tangential to the surface and $x$ and $r$ for the axial and radial cylindrical coordinates. Additional constraints include the balance of stresses on the surface along the normal $\boldsymbol{n}$ and tangential $\boldsymbol{t}$ directions

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle p=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{n}-\frac{2}{\unicode[STIX]{x03C0}Re}\boldsymbol{n}\boldsymbol{\cdot }\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D707}}\boldsymbol{\cdot }\boldsymbol{n}-\frac{1}{2\unicode[STIX]{x03C0}^{2}\unicode[STIX]{x1D6F1}^{1/2}}[E_{n}^{o2}-\unicode[STIX]{x1D700}E_{n}^{i2}+(\unicode[STIX]{x1D700}-1)E_{t}^{2}], & \displaystyle\end{eqnarray}$$

(3.5) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{t}\boldsymbol{\cdot }\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D707}}\boldsymbol{\cdot }\boldsymbol{n}=\frac{Re}{2}\unicode[STIX]{x1D70E}E_{t}, & \displaystyle\end{eqnarray}$$

( $\unicode[STIX]{x1D749}_{\unicode[STIX]{x1D707}}$ is the viscous stress tensor), the kinematic velocity condition and the jump of the normal component of the electric field at the surface

(3.6) $$\begin{eqnarray}\displaystyle \boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{n}=0,\quad \unicode[STIX]{x1D70E}=\frac{2}{\unicode[STIX]{x03C0}\unicode[STIX]{x1D6F1}^{1/2}}(E_{n}^{o}-\unicode[STIX]{x1D700}E_{n}^{i}). & & \displaystyle\end{eqnarray}$$

A velocity field is prescribed at the inlet boundary (we use a homogeneous normal velocity component injecting a flow rate $Q$ ), while the Taylor potential is imposed as the electrostatic far-field boundary condition (Taylor Reference Taylor1964; Gañán-Calvo Reference Gañán-Calvo1997; Higuera Reference Higuera2003):

(3.7) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D719}_{T}(x,r)=-4.23\unicode[STIX]{x1D6F1}^{1/4}(x^{2}+r^{2})^{1/4}P_{1/2}\left(\frac{x}{\sqrt{x^{2}+r^{2}}}\right). & & \displaystyle\end{eqnarray}$$

Here, $P_{1/2}$ is the Legendre function of the first kind of degree $\frac{1}{2}$ . The electrostatic boundary condition is imposed on a circumference centred at the apparent vertex of the meniscus, with a dimensionless radius of $200\cdot \unicode[STIX]{x1D700}^{1/6}$ . The operational conditions simulated include the dielectric constants of tributyl phosphate and propylene carbonate, 8.91 and 64.9; two Reynolds numbers for each dielectric constant (0.290 and 1.12 for $\unicode[STIX]{x1D700}=8.91$ , and 0.396 and 0.740 for $\unicode[STIX]{x1D700}=64.9$ ); and a wide range of flow rates for each Reynolds number, down to the minimum observed in experiments.

Figure 4 compares the current obtained from the model with experimental data for several solutions in table 1. The total current $\tilde{I}$ is the sum of the bulk conduction current $\tilde{I} _{B}(\tilde{x})=\int _{0}^{\tilde{R}}K\tilde{E}_{x}^{i}2\unicode[STIX]{x03C0}\tilde{r}\,\text{d}\tilde{r}$ , and the surface current $\tilde{I} _{S}(\tilde{x})=2\unicode[STIX]{x03C0}\tilde{R}\tilde{\unicode[STIX]{x1D70E}}\tilde{v}_{s}$ . The total current is constant along the axial coordinate due to conservation of charge. We use the scale $I_{o}=(\unicode[STIX]{x1D700}_{o}\unicode[STIX]{x1D6FE}^{2}/\unicode[STIX]{x1D70C})^{1/2}$ to illustrate the scaling law $\tilde{I} /I_{o}\propto \unicode[STIX]{x1D6F1}^{1/2}$ , which can be justified as an asymptote for large $\unicode[STIX]{x1D6F1}$ (Fernández de la Mora & Loscertales Reference Fernández de la Mora and Loscertales1994; Gañán-Calvo et al. Reference Gañán-Calvo, Dávila and Barrero1997; Higuera Reference Higuera2003). The numerical results and the experimental data compare well for both liquids, supporting the validity of the model. A similar comparison by Gamero-Castaño & Magnani (Reference Gamero-Castaño and Magnani2019) yielded a good but inferior match, probably because in this case the electrosprays were operated under vacuum (Gamero-Castaño Reference Gamero-Castaño2010), with flow rates indirectly measured via time of flight. This technique is affected by uncertainties associated with the spreading of the beam and losses of mechanical energy in the cone jet (Grustan-Gutierrez & Gamero-Castaño Reference Grustan-Gutierrez and Gamero-Castaño2017), and is less reliable than the direct flow rate measurement employed in this work. The numerical solution is well fitted by the linear functions $\tilde{I} /I_{o}=2.66\unicode[STIX]{x1D6F1}^{1/2}-0.48$ for tributyl phosphate, and $\tilde{I} /I_{o}=2.63\unicode[STIX]{x1D6F1}^{1/2}-5.46$ for propylene carbonate. The Reynolds number is the weakest contributor in the $I(\unicode[STIX]{x1D6F1},Re,\unicode[STIX]{x1D700})$ relationship, to the point that the current of a cone jet is generally regarded as independent of the viscosity of the liquid. On the other hand the current does depend on the dielectric constant, e.g. in the $\tilde{I} /I_{o}\propto \unicode[STIX]{x1D6F1}^{1/2}$ scaling law this dependency appears in the form of a negative $y$ -intercept that becomes more significant at increasing dielectric constant.

Figure 4.

Current versus flow rate, comparison between numerical solution and experiments. The data follow the scaling law $\tilde{I} /I_{o}\sim \unicode[STIX]{x1D6F1}^{1/2}$ , with $I_{o}=(\unicode[STIX]{x1D700}_{o}\unicode[STIX]{x1D6FE}^{2}/\unicode[STIX]{x1D70C})^{1/2}$ . The linear fits are for the numerical data.

4 Analysis of the numerical solution and discussion

We next analyse mechanisms that may cause the minimum flow rate, starting with the hypothesis that forces impeding the flow become dominant at decreasing flow rate. We will show that the algorithm is able to find numerical solutions at flow rates below the experimental minimum, especially in the case of the less polar liquid tributyl phosphate. Therefore the evidence for positive identification of the mechanism will be indirect and the theoretical conclusions somewhat speculative. The order of the analysis is reduced by using the following differential balance between jet kinetic power, pressure flow work, viscous stress flow work and the linear powers injected by the tangential electric stress on the surface and removed by viscous dissipation from the bulk (see appendix A for its derivation),

(4.1) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x1D70C}Q}{2}\frac{\text{d}\langle \tilde{v}^{2}\rangle }{\text{d}\tilde{x}} & = & \displaystyle \sqrt{1+\tilde{R}^{\prime 2}}\tilde{I} _{S}\tilde{E}_{t}-Q\frac{\text{d}\langle \tilde{p}\rangle }{\text{d}\,\tilde{x}}\nonumber\\ \displaystyle & & \displaystyle +\,\frac{\text{d}}{\text{d}\,\tilde{x}}\int _{0}^{\tilde{R}}(\tilde{\unicode[STIX]{x1D70F}}_{\unicode[STIX]{x1D707}}\boldsymbol{\cdot }\boldsymbol{i}_{x})\boldsymbol{\cdot }\tilde{\boldsymbol{v}}2\unicode[STIX]{x03C0}\tilde{r}\,\text{d}\tilde{r}-\int _{0}^{\tilde{R}}(\tilde{\unicode[STIX]{x1D70F}}_{\unicode[STIX]{x1D707}}:\unicode[STIX]{x1D735}\tilde{\boldsymbol{v}})2\unicode[STIX]{x03C0}\tilde{r}\,\text{d}\tilde{r}.\end{eqnarray}$$

The averages stand for variables integrated in the cross-section of the cone jet and weighted with the axial component of the velocity, e.g. $\langle \tilde{p}\rangle =(1/Q)\int _{0}^{\tilde{R}}\tilde{p}\tilde{v}_{x}2\unicode[STIX]{x03C0}\tilde{r}\,\text{d}\tilde{r}$ . The dimensionless form of equation (4.1) is

(4.2) $$\begin{eqnarray}\displaystyle \underbrace{\frac{\unicode[STIX]{x1D6F1}_{Q}^{1/2}}{2\unicode[STIX]{x03C0}^{2}}\frac{\text{d}\langle v^{2}\rangle }{\text{d}x}}_{T1} & = & \displaystyle \underbrace{\sqrt{1+R^{\prime 2}}I_{s}E_{t}}_{T2}\underbrace{-\frac{\text{d}\langle p\rangle }{\text{d}x}}_{T3}\nonumber\\ \displaystyle & & \displaystyle \underbrace{+\frac{1}{\unicode[STIX]{x03C0}Re}\frac{\text{d}}{\text{d}x}\int _{0}^{R}(\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D707}}\boldsymbol{\cdot }\boldsymbol{i}_{x})\boldsymbol{\cdot }\boldsymbol{v}2\unicode[STIX]{x03C0}r\,\text{d}r}_{T4}\underbrace{-\frac{1}{\unicode[STIX]{x03C0}Re}\int _{0}^{R}(\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D707}}\boldsymbol{ : }\unicode[STIX]{x1D735}\boldsymbol{v})2\unicode[STIX]{x03C0}r\,\text{d}r}_{T5}.\end{eqnarray}$$

Strictly a power balance, (4.2) can also be interpreted as a balance of the average of forces projected on the direction of the flow. As the cone jet becomes slender (4.2) reduces to the one-dimensional momentum equation derived by Gañán-Calvo (Reference Gañán-Calvo1999) and used by Ponce-Torres et al. (Reference Ponce-Torres, Rebollo-Muñóz, Herrada, Gañán-Calvo and Montanero2018) to analyse their numerical solutions.

Figure 5.

Numerical solution for $\unicode[STIX]{x1D700}=8.91$ (tributyl phosphate), $Re=0.29$ and two flow rates, $\unicode[STIX]{x1D6F1}=5.44,31.6$ . The top charts plot the position of the surface, its second derivative, the conduction and surface currents and the ratio between the surface current and its value in equilibrium. The bottom charts plot the components of the balance (4.2).

Figure 5 illustrates the numerical solution for tributyl phosphate, a Reynolds number of 0.29 and two flow rates, $\unicode[STIX]{x1D6F1}=5.44$ and $31.6$ . The lower value is slightly smaller than the experimental minimum at this Reynolds number, $\unicode[STIX]{x1D6F1}_{min}=7.30$ . The top charts show the position of the surface, its second derivative $R^{\prime \prime }(x)$ , the bulk conduction current, the surface current and the ratio between the surface current and its value if the surface charge were in equilibrium, $I_{S}/I_{S,eq}=1-\unicode[STIX]{x1D700}E_{n}^{i}/E_{n}^{o}$ . In equilibrium no charge is injected on the surface from the bulk, i.e. $E_{n}^{i}=0$ and $\unicode[STIX]{x1D70E}_{eq}=2E_{n}^{o}/(\unicode[STIX]{x03C0}\unicode[STIX]{x1D6F1}^{1/2})$ . The bottom charts plot the components of the balance (4.2). The conical surface turns into a jet in an interval centred at the maximum of $R^{\prime \prime }$ , and using its position $x_{o}$ as origin of the axial coordinate helps comparing the solutions for different $\{\unicode[STIX]{x1D6F1},Re,\unicode[STIX]{x1D700}\}$ states. In the cone, as the liquid moves downstream the increasing surface velocity reduces the surface charge below its equilibrium value, and an electric field $E_{n}^{i}$ develops to inject charge on the surface and reduce the imbalance. Consequently $I_{B}$ decreases and $I_{S}$ increases in equal amounts (the total current $I_{B}(x)+I_{S}(x)$ remains constant). The two currents become equal at the cross-over point $x_{c}$ , located a few units downstream of $x_{o}$ . The interval within which the conduction and surface currents significantly reverse their values defines the transition region of the cone jet. In addition we define the shorter geometry transition as the interval where the cone evolves into a jet; identifiable by the significant value of $R^{\prime \prime }$ , it occupies the initial section of the transition region. The transition region is characterized by the absence of surface charge equilibrium (see $I_{S}/I_{S,eq}$ curves), to a degree that depends on the ratio between the electric relaxation time and the fluid residence time. This ratio can be cast in the form $\unicode[STIX]{x1D700}/\sqrt{\unicode[STIX]{x1D6F1}}$ (Gamero-Castaño & Magnani Reference Gamero-Castaño and Magnani2019), which naturally arises when rewriting equation (3.2) in terms of the surface current and the equilibrium surface current (Fernández de la Mora & Loscertales Reference Fernández de la Mora and Loscertales1994)

(4.3) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x1D700}}{\unicode[STIX]{x03C0}\unicode[STIX]{x1D6F1}^{1/2}}\frac{\text{d}I_{S}}{\text{d}t}=I_{S,eq}-I_{S}, & & \displaystyle\end{eqnarray}$$

with $\text{d}t=\sqrt{1+R^{\prime 2}}(\text{d}x/v_{s})$ . The departure from equilibrium increases with $\unicode[STIX]{x1D700}/\sqrt{\unicode[STIX]{x1D6F1}}$ , so that for a liquid with low dielectric constant like tributyl phosphate the departure from full relaxation is small for most flow rates realized experimentally, and relatively significant near the minimum flow rate. The absence of charge relaxation at flow rates that are stable in experiments invalidates the hypothesis behind the criterion (1.2) for the minimum flow rate. The bottom charts in figure 5 show that the tangential electric stress and viscous stress are the dominant forces accelerating and impeding the flow respectively. The pressure increases monotonically in the direction of the flow, also impeding the motion for both flow rates. Only at very high flow rates, e.g. $\unicode[STIX]{x1D6F1}\gtrsim 66.5$ for $Re=0.290$ , the pressure develops a negative peak in the geometry transition, followed by a gradual increase to positive values. This is due to the increasing importance of the electrostatic pull, i.e. the $-E_{n}^{o2}/(2\unicode[STIX]{x03C0}^{2}\unicode[STIX]{x1D6F1}^{1/2})$ term in equation (3.4), and the scaling $E_{n}^{o}\sim \unicode[STIX]{x1D6F1}^{1/2}$ near the current cross-over point (Gañán-Calvo Reference Gañán-Calvo2004; Gamero-Castaño & Magnani Reference Gamero-Castaño and Magnani2019). The derivative of the jet kinetic power, and therefore the fluid’s acceleration, is positive throughout the transition region and increases monotonically downstream for both flow rates. Thus the average acceleration profile does not provide an indication that the flow is decelerating, much less stopping, which could make the flow unstable.

Figure 6.

Linear electric power accelerating the flow ( $T2$ ) and total linear viscous power impeding the flow ( $T4+T5$ ) as a function of the flow rate, computed at both the current cross-over point $x_{c}$ and the geometry transition centre $x_{o}$ . Calculations are for $\unicode[STIX]{x1D700}=8.91$ (tributyl phosphate), $Re=0.29$ and $Re=1.12$ .

The acceleration of the flow increases monotonically downstream for all simulated flow rates of tributyl phosphate, including those below the minimum, and stays positive throughout the transition region. However the acceleration profiles do become lower at decreasing flow rate. Since a sufficiently low value of the acceleration, rather than zero, may be triggering the instability responsible for the minimum flow rate, a balance between opposing stresses could still yield the experimental dependency $\unicode[STIX]{x1D6F1}_{min}(Re)$ observed in figure 3. Figure 6 tests this hypothesis by plotting the linear electric power $T2$ accelerating the fluid and the net linear viscous power impeding it, $T4+T5$ , at two representative axial locations, the maximum of $R^{\prime \prime }$ and the current cross-over point, i.e. in the initial section of the transition region and at its centre. The data are for tributyl phosphate and two Reynolds numbers, 0.290 and 1.12. The impediment to the flow by viscous stresses is larger at $x_{o}$ than at $x_{c}$ , and remains approximately constant with flow rate. On the other hand the accelerating electric power at $x_{c}$ scales as $\unicode[STIX]{x1D6F1}_{Q}^{0.22}$ for $Re=0.290$ and $\unicode[STIX]{x1D6F1}_{Q}^{0.34}$ for $Re=1.12$ , while it is smaller and remains relatively constant at $x_{o}$ . Thus, a balance of the two dominant terms at $x_{o}$ suggests that it has no effect on the stability of the flow, while the balance at $x_{c}$ suggests the law $\unicode[STIX]{x1D6F1}_{min}\propto Re^{-a}$ with $a$ between 2.9 and 4.5. Both scenarios differ significantly from the experimental law $\unicode[STIX]{x1D6F1}_{min}\propto Re^{-1.1}$ .

Figure 7.

Streamlines (with values of the streamfunction) and contour map of the norm of the velocity for $\unicode[STIX]{x1D700}=8.91$ (tributyl phosphate), $Re=0.29$ and two flow rates, $\unicode[STIX]{x1D6F1}=5.44$ and $\unicode[STIX]{x1D6F1}=31.6$ .

The previous analysis suggests that the minimum flow rate is not determined by a near balance between opposing stresses. Because the integration of the variables leading to (4.2) may be suppressing two-dimensional features associated with flow deceleration and instability, figure 7 provides details of the flow field including streamlines (the streamfunction varies from 0 at the axis to $1/2$ at the free surface), and a contour map of the norm of the velocity. The plots extend upstream to areas where the terms in (4.2) are negligible compared to their maximum values. The streamlines are similar for both flow rates and nearly parallel to the free surface, while the norm of the velocity is almost constant at fixed $x$ . Thus, the flow field in the transition region is quasi-unidirectional (parallel to the free surface), changes little with the flow rate and there is no local flow deceleration nor indication that the flow may become unstable. Note also that the constant velocity norm at fixed $x$ everywhere in the transition region, including areas where the geometry is not slender, increases the accuracy of one-dimensional jet models like that of Gañán-Calvo (Reference Gañán-Calvo1999).

Figure 8.

Numerical solution for $\unicode[STIX]{x1D700}=64.9$ (propylene carbonate), $Re=0.74$ and two flow rates, $\unicode[STIX]{x1D6F1}=30.7,200.8$ . The top charts plot the position of the surface, its second derivative, the conduction and surface currents and the ratio between the surface current and its value in equilibrium. The bottom charts plot the components of the balance (4.2).

Figure 8 shows the main features of the transition region for propylene carbonate, a Reynolds number of 0.740 and dimensionless flow rates 30.7 and 200.8. The lower flow rate is near the minimum experimental one, $\unicode[STIX]{x1D6F1}_{min}=19.1$ . The much larger dielectric constant of propylene carbonate has important consequences in the transition region. First, because of the larger ratio between the electric relaxation and the fluid residence times it takes a longer distance for the surface charge to equilibrate, making the transition region longer; and the surface charge falls further from equilibrium (compare the $I_{S}/I_{S,eq}$ profiles of both liquids), making it possible for the conical shape to extend farther towards the axis before turning into a jet. Consequently the profile of the geometry transition is much sharper, i.e. the second derivative $R^{\prime \prime }$ exhibits a higher peak and is narrower for large values of $\unicode[STIX]{x1D700}/\sqrt{\unicode[STIX]{x1D6F1}}$ (see $R^{\prime \prime }$ profile for $\unicode[STIX]{x1D6F1}=30.7$ ). Second, at low flow rates the pressure gradient in (4.2) becomes the dominant accelerating term in the geometry transition, and rapidly switches to an impeding role near the current cross-over point. This is due to the $-(\unicode[STIX]{x1D700}-1)E_{t}^{2}/(2\unicode[STIX]{x03C0}^{2}\unicode[STIX]{x1D6F1}^{1/2})$ component of the electrostatic stress in equation (3.4) and the scaling $E_{t}\sim \text{constant}$ near the current cross-over point. Both combined cause a large negative pressure peak near the current cross-over at high values of $\unicode[STIX]{x1D700}/\sqrt{\unicode[STIX]{x1D6F1}}$ . The viscous stress remains the dominant impeding force in the geometry transition, while the electrostatic tangential stress becomes the dominant accelerating force downstream. Similarly to tributyl phosphate the acceleration of the fluid is always positive in the transition region, although at low flow rates the profile is not monotonic due to the rapid sign reversal of the pressure gradient. A balance of dominant forces at the current cross-over point cannot explain the dependence of the minimum flow rate on viscosity, because the dominant balance is between the electric tangential stress and the pressure gradient and does not include the viscous stress. Such balance could occur upstream in the geometry transition region where the viscous stress is the dominant impeding force, but in this area the accelerating pressure gradient is most intense, the acceleration of the flow exhibits a local maximum, and it is unlikely that the instability develops under the most favourable flow condition. In summary, the numerical solution suggests that the experimental dependence of the minimum flow rate on viscosity is not associated with a balance of stresses.

Figure 9.

Power dissipated ( $P_{\unicode[STIX]{x1D707}}$ ) and power injected by the tangential electric stress ( $P_{E}$ ) up to the current cross-over point, as a function of the flow rate. (a) Propylene carbonate, $Re=0.396$ and $Re=0.74$ ; (b) tributyl phosphate, $Re=0.29$ and $Re=1.12$ .

The dependence of the minimum flow rate on the Reynolds number requires an alternative instability mechanism with the participation of the viscous stress. The profiles in figures 5 and 8 hint at one possibility: viscous dissipation, which predominantly takes place in the geometry transition, first becomes significant upstream of the work done by the electric stress. Since the tangential electric stress is the only external force doing work on the system, it must couple with the pressure and velocity fields so that the energy that needs to be dissipated can be transferred upstream. Note also that the viscous dissipation is directly related to the geometry of the cone jet: for this quasi-one-directional flow $\unicode[STIX]{x2202}u/\unicode[STIX]{x2202}x$ is the larger component of the stress tensor, $u\approx 1/R^{2}$ , and $T5\approx -(1/Re)(8R^{\prime 2}/R^{4})$ . Figure 9 explores the relationship between the minimum flow rate and the injection and dissipation of mechanical energy. It shows the power dissipated up to the current cross-over point, $P_{\unicode[STIX]{x1D707}}(x_{c})=-\int _{-\infty }^{x_{c}}T5\,\text{d}x$ , and the power injected by the tangential electric stress, $P_{E}(x_{c})=\int _{-\infty }^{x_{c}}T2\,\text{d}x$ , as functions of the dimensionless flow rate. Figure 9(a) is for tributyl phosphate and $Re=0.290,1.12$ , and figure 9(b) is for propylene carbonate and $Re=0.396,0.740$ . The location of the current cross-over is significant for this analysis because it approximately bounds the region where the tangential electric stress, the viscous stress and the pressure gradient in the case of large $\unicode[STIX]{x1D700}/\sqrt{\unicode[STIX]{x1D6F1}}$ are comparable. Farther downstream $R^{\prime }$ decays rapidly decreasing the viscous dissipation density, while the bulk conduction current becomes small compared to the surface current. Both trends increase the conversion of the electric work into kinetic energy. At high flow rates the power injected up to $x_{c}$ , i.e. within the geometry transition, is much larger than the power dissipated. However as the flow rate decreases the difference becomes smaller and the two functions intersect at a flow rate that is comparable to the minimum one: e.g. the tributyl phosphate curves intersect at $\unicode[STIX]{x1D6F1}=20.1$ for $Re=0.290$ and $\unicode[STIX]{x1D6F1}=3.1$ for $Re=1.12$ , while the minimum flow rates are 7.3 and 1.7 respectively; and the intersection for propylene carbonate is at $\unicode[STIX]{x1D6F1}=31.0$ for $Re=0.396$ and $\unicode[STIX]{x1D6F1}=17.0$ for $Re=0.740$ , while the minimum flow rates are 26.2 and 19.1 respectively. Thus, figure 9 shows that as the flow rate decreases a larger fraction of the energy that is dissipated is being injected into the system farther downstream. The choice of $x_{c}$ is somewhat arbitrary, and the fitting parameters of the curves in figure 9 depend on this choice. Thus, we do not expect these fitting parameters to reproduce accurately the experimental dependency between the minimum flow rate and the Reynolds number. Furthermore a fraction of the work done by the electric stress is converted into kinetic energy. Despite this, the intersections of $P_{E}(x_{c})$ and $P_{\unicode[STIX]{x1D707}}(x_{c})$ yield $\unicode[STIX]{x1D6F1}_{min}=3.6Re^{-1.4}$ for tributyl phosphate and $\unicode[STIX]{x1D6F1}_{min}=12.7Re^{-0.96}$ for propylene carbonate, which are not very different from the experimental relationships in figure 3 and reproduce the trend of a minimum flow rate increasing with dielectric constant. The separation between the regions where power is injected and dissipated could be the mechanism making the flow unstable, as a perturbation of the flow between these two regions will alter the transfer of energy that must be dissipated, and which is fixed by the position of the surface in the geometry transition. Furthermore, the instability of the system is likely to increase with both the amount of energy that needs to be transferred and with the transfer distance, i.e. the instability is likely to increase at decreasing flow rate.

Gamero-Castaño (Reference Gamero-Castaño2019) has analysed both the ohmic and viscous dissipation taking place in cone jets to find the conditions under which self-heating is important and isothermal operation breaks down. The scaling law $P_{\unicode[STIX]{x1D6FA}}\cong 7.30\unicode[STIX]{x1D6F1}^{0.41}$ fits well the total ohmic dissipation in the transition region for both tributyl phosphate and propylene carbonate. The total viscous dissipation exhibits a weaker dependence on the minimum flow rate (see figure 9), and is linear in $Re^{-1}$ : $P_{\unicode[STIX]{x1D707}}\cong 2.47Re^{-1}+2.61$ for tributyl phosphate at $\unicode[STIX]{x1D6F1}=31.6$ , and $P_{\unicode[STIX]{x1D707}}\cong 1.84Re^{-1}+1.77$ for propylene carbonate at $\unicode[STIX]{x1D6F1}=101$ . Although the ohmic dissipation exceeds the viscous dissipation under the conditions investigated in this article, it has no effect in the balance of mechanical energy used in figure 9 because none of the mechanical work done by the electric stress is dissipated ohmically. Instead, ohmic dissipation is directly associated with the energy injected by the electric field to drive conduction current in the bulk of the liquid, none of which is transferred in the form of mechanical work (see appendix A).

We have adopted throughout the article the point of view that the minimum flow rate is triggered by a hypothetical flow instability. This is supported by the work of Ponce-Torres et al. (Reference Ponce-Torres, Rebollo-Muñóz, Herrada, Gañán-Calvo and Montanero2018) who obtain numerical solutions of the leaky dielectric model at flow rates below the minimum one, and which upon further investigation are found to be unstable. A similar inference results from our simulations of tributyl phosphate, for which we find solutions at flow rates below the minimum experimental one. On the other hand Higuera (Reference Higuera2017) argues that the minimum flow rate is an intrinsic value that the electric shear stress can pull under important charge relaxation effects, and that no solution may exist at lower flow rates. In this regard we note that convergence in our algorithm becomes increasingly difficult at decreasing flow rate near the minimum one, especially for the more polar liquid propylene carbonate. For example we have been unable to obtain solutions down to the minimum flow rate for propylene carbonate at $Re=0.74$ . However convergence depends on the initial shape of the profile, and it is not possible for us to say whether solutions cease to exist below a critical flow rate.