In this work, non-reactive molecular dynamic simulations were conducted to determine the surface tension of water as a function of the concentration of the dissolved gaseous molecules (O2), which would in turn help to predict the pressure inside the nanobubbles under supersaturation conditions. Knowing the bubble pressure is a prerequisite for understanding the mechanisms behind the spontaneous combustion of the H2/O2 gases inside the nanobubbles. First, the surface tension of pure water was determined using the planar interface method and the Irving and Kirkwood formula. Next, the surface tension of water containing four different supersaturation concentrations (S) of O2 gas molecules was computed considering the curved interface of a nanobubble. The surface tension of water was found to decrease with an increase in the supersaturation ratio or the concentration of the dissolved O2 gas molecules.

Combustion-based small-scale energy devices have applications in many fields, e.g., portable power devices, sensors, actuators, unmanned aerial microvehicles, microthrusters, and micro-heating devices.1–3 At smaller scales, the surface-area-to volume ratio of the combustor increases, and the characteristic inertia time decreases. This leads to reduced mixing, increased heat loss from the wall, and strong wall-flame kinetic interactions.4 The successful development of small-scale power generation systems continues to face significant challenges. One of the biggest difficulties is that as a result of the increased surface-area-to-volume ratio at smaller scales, combustion is found difficult to be initiated and sustained because heat loss from the wall exceeds heat release from the chemical reactions.5 Loyal to the scaling law, we would expect combustion to be impossible at nanoscale because heat loss would profoundly dominate chemical reactions. However, this is not always the case. Svetovoy et al.6,7 accidentally discovered that the stoichiometric H2/O2 mixture in nanobubbles could be ignited spontaneously with bubble sizes smaller than 150 nm. These nanobubbles were produced from short-time (1-100 μs) water electrolysis by applying high-frequency alternating sign voltage pulses, which resulted in the H2 and the O2 gas production above the same electrode. Motivated by this, Jain et al.8 developed an independent experiment, similar to that of Svetovoy et al.,6,7 but with more quantitative determination of the temperature change of the electrode surface, which was achieved by using an in-built micro-fabricated Pt resistance thermal sensor. They found that the combustion only occurred for frequencies greater than 15 kHz and that the amount of heat released was maximized for a duty cycle of 0.5, corresponding to the stoichiometric H2/O2 gas production. These were consistent with the findings of Svetovoy et al.,6,7 suggesting that the combustion had indeed taken place inside the nanobubbles. Recently, in a review paper by Svetovoy et al.,9 the various experimental work conducted by their group in relation to the combustion inside the nano/micro-bubbles was summarized. Moreover, the authors proposed the combustion process to be surface-dominated and to occur at room temperature. Thus, knowing the bubble pressure is a significant first step before the reaction kinetics could be understood.

As a result of the capillary effect, the pressure inside a nanobubble is much higher as compared to the pressure inside a micro-sized bubble, which implies that the gases inside the nanobubble cannot be in equilibrium with the atmosphere, since larger pressure means larger chemical potential, and consequently the nanobubble should dissolve (in a saturated environment).10 This contradicts the longevity and the apparent stability of nanobubbles observed experimentally.11–36 The nanobubbles can be divided into two categories: surface nanobubbles and bulk nanobubbles. The surface nanobubbles are gas-filled cavities on a surface in the form of a spherical cap with the uniform radius of curvature and geometrical center located beneath the surface, whereas the bulk nanobubbles are spherical shaped gas-filled cavities located in the interior of a liquid. The first experimental evidence of bulk nanobubbles was reported by Johnson et al.,11 where the bubbles were produced by shear in seawater and were observed to be stable for more than 22 hrs. Following this work, a number of other studies were conducted to understand the generation and stabilization of these bulk nanobubbles.12–22 As far as the surface nanobubbles are considered, Parker et al.23 first reported their existence to explain the long-range attraction between the hydrophobic surfaces. After which, many AFM (atomic force microscopy) measurements have confirmed the existence of these spherical-caped gaseous nanobubbles at solid-liquid interfaces.24–36 Refer to the reviews given by Alheshibri et al.37 and Lohse38 for a more detailed overview on the bulk and the surface nanobubbles. In our previous experimental work5 and in studies performed by Svetovoy et al6,7,9 both kind of nanobubbles (surface and bulk) were produced during the water electrolysis processes.

Three hypotheses can explain the contradiction of the experimental evidence of nanobubble existence and the paradoxical L-Y (Laplace-Young) equation which suggests that nanobubbles cannot be stable: (1) The L-Y equation is not valid at nanoscale; (2) The nanobubbles are in equilibrium but in a supersaturated solution; (3) “The traffic jam effect.”39 

  1. The L-Y equation, which assumes the application of macroscopic, continuum thermodynamics on the nanoscopic length scales, is derived from balancing the mechanical forces on the bubble surface. Thus, if the measured stability of nanobubbles is taken as a proof to invalidate the L-Y equation, then the applications of similar continuum theories on the nanoscopic length scales such as, the application of macroscopic hydrodynamics to microfluidic flow and the DLVO (Derjaguin, Landau, Verwey, and Overbeek) theory, need to be reconsidered. However, there are a number of experimental and computational work done to support the quantitative application of such continuum theories on the nanoscopic scales with the nanobubbles being the only exception.10 Thus the L-Y equation, which is nothing but a minimization of the total free energy with respect to the bubble radius can be assumed to be valid even for the nanobubbles.

  2. For a nanobubble to be completely stable in a liquid, it needs to be in both mechanical and chemical equilibrium simultaneously. The application of the L-Y equation assumes only mechanical equilibrium, thus for the bubble to be in chemical equilibrium, the chemical potential of the gases inside and outside the bubble must also be equal i.e. the liquid must be super-saturated with the bubble gases. The supersaturation (S) is defined as the ratio of the concentration of the gas (C) in the liquid to the saturated concentration of the gas, at 1 atm, (C*) i.e. S = C/C*. Moody40 and He41 conducted monte-carlo simulations to calculate the surface tension of supersaturated Lennard-Jones planar liquid-vapor interfaces. They found the surface tension value to decrease with increasing supersaturation and to vanish at the vapor spinodal. Thus, the surface tension of water needs to be modified according to the local supersaturation value in order to get the correct pressure difference across the nanobubble. All of the previous studies40,41 relating to supersaturations have been performed for pure systems. However, in this work heterogeneous systems were considered where molecular dynamic simulations were conducted to determine the surface tension of water as a function of the concentration of the dissolved gaseous molecules under supersaturation conditions. Here only O2 molecules were considered but the qualitative behavior of the surface tension value with other types of gases like H2 will be similar. Thus, the main purpose was to determine the internal pressure of the nanobubbles within which spontaneous combustion of the H2/O2 gases takes place.

  3. A nanobubble can also be stable in a closed container where the gases are not allowed to leave the liquid into the atmosphere. Due to the initial diffusion of the gases out of the nanobubbles, the liquid surrounding the nanobubbles could become supersaturated and thus an equilibrium could be established that stabilizes the nanobubbles. However, if the container is opened, then the gases leave the liquid, allowing the nanobubbles to completely dissolve. This effect is called “the traffic jam effect”, as studied by Weijs.39 However, considering our experimental conditions,5 the explanation (2) is more relevant as the liquid droplet was exposed to the atmosphere with the liquid being continuously supersaturated with the gases produced from the water electrolysis.

The simulations were performed using LAMMPS,42 an open source MD simulation code developed by Sandia National Lab. The simulation system consists of three types of interactions: (1) H2O–H2O Interactions: The water molecules were modeled using the extended simple point charge (SPC/E) force field,43 which considers both the columbic and the dispersive Lennard-Jones (LJ) interactions. A partial charge of −0.8476|e| and 0.4238|e| was used for the oxygen and the hydrogen atom respectively, where |e| refers to the electron charge unit. A rigid tetrahedral angle (H–O–H) of 109.47o and a bond length (O-H) of 0.1 nm was assumed. In addition, the SHAKE algorithm44,45 was used to rigidify the water molecules; (2) O2–O2Interactions: The oxygen molecules were modeled as linear with the partial charges on the oxygen atoms being assumed to be zero. The equilibrium bond length (O-O) was set to 0.121 nm. Both the L-J and the columbic interactions were considered similar to the H2O-H2O interactions; (3) H2O–O2Interactions: The interactions between the water and the oxygen molecules were obtained again by adding the columbic and the L-J contributions. The values of σH2O-O2 and εH2O-O2 terms were calculated by using the Lorentz-Berthelot mixing rules:46,47

(1)

The L-J parameters for each atom type are given in Table I:48 

TABLE I.

LJ parameters.

ParameterH2O (O atom)H2O (H atom)O2 (O atom)
σ (A°) 3.166 3.094 
ε (kcal/mol) 0.1553 0.09538 
ParameterH2O (O atom)H2O (H atom)O2 (O atom)
σ (A°) 3.166 3.094 
ε (kcal/mol) 0.1553 0.09538 

In order to validate the SPC/E water model and the simulation parameters, the surface tension of pure water was first determined. A planar interface was used to calculate the surface tension according to the Irving and Kirkwood formula.49 The simulation box, as shown in Fig. 1, was divided into three slabs of sizes (7 x 7 x 7 nm) with the middle slab consisting of N = 12,000 water molecules based on the water density at P = 1 atm and T = 300 K. Periodic boundary conditions were used in all the three directions and a time step of 2 fs was used. After the initial set-up of the computational domain, the system equilibration to 300 K and 1 atm using the Nose-Hoover thermostat50,51 and barostat52 with a relaxation time of 200 fs and 1000 fs, respectively, was performed for 2 ns. Then the production run was performed under the NVT (constant number, volume and temperature) ensemble for another 2 ns. For all the simulations, a high L-J cut-off value of 2.5 nm was used to ensure that errors introduced into the surface tension calculations due to the L-J tail-corrections were minimum. A surface tension value of 59.83 ± 2 mN/m was obtained similar to the value of 63.6 mN/m obtained by Vega et al.53 using the SPC/E rigid water model. In a study performed by Vega et al.,53 much better agreement with the experimental surface tension value of 71.3 mN/m was obtained by the SPC/E, TIP4P/Ew, and TIP4P/2005 models than the SPC, TIP3P, and TIP4P water models.53 Thus, because of the wide use and simplicity, all of the following surface tensions were conducted using the SPC/E rigid water model.

FIG. 1.

A snapshot of the simulation box.

FIG. 1.

A snapshot of the simulation box.

Close modal

In the nanobubble experimental study, steady faradic currents up to 0.3 mA were obtained for a 20 μm x 20 μm sized electrode surface.8 Using the diffusion length and the faradic current equation, the amount of super-saturation in the liquid water surrounding the electrodes could be estimated. The diffusion length (lO2) is given by the following equation:7 

(2)

where DO2 is the diffusion coefficient of the oxygen molecule6 (2x10−9 m2/s) and t is half the time period of a complete cycle. The concentration of the dissolved oxygen molecules (CO2) can be calculated using the faradic current equation7 as follows:

(3)

where I is the faradic steady state current, e is the electron charge value and A is the cross-sectional area of the electrode surface. Using the steady state current value of 0.3 mA, 20 μmx20 μm electrode surface and t = 33 μs (for a 15 kHz alternating polarity pulse), the CO2 value comes out to be 1.5x1020 cm−3. The relative supersaturation (S), which is the ratio of the gas concentration (CO2) to the saturated concentration, CO2* = 7.7x1017 cm−3, (at P = 1 atm and T = 300 K)7 comes out to be around 200. Thus, in the present MD simulations, the water system was supersaturated with O2 gas molecules with S values ranging from 1-850, which covers most of the saturation levels achieved in the experiment.

Initially, the O2 and the H2O molecules were randomly added to the simulation box (7 nm x 7 nm x 7nm). The number of initial O2 gas molecules added was based on the supersaturation value desired. First the system equilibration to 300 K and 1 atm was performed for 2 ns. The relaxation time for the Nose-Hoover thermostat and barostat was kept same as that of the pure water simulation i.e. 200 fs and 1000 fs, respectively. A time step of 2 fs was used with the 3d bin sizes being 7 nm x 7 nm x 5 A°. The L-J-cut off value was 2.5 nm and periodic boundary conditions were used in all the three directions. After the equilibration, a cavity of size 2.5 nm was created at the center of the simulation box (by deleting some of the water molecules). The simulation was then run under the NVT conditions (300 K) until the bubble radius attains a steady state value, after which the pressure inside and outside of the bubble was calculated. The surface tension value was calculated by assuming the validity of the Laplace-Young equation and by using the ΔP across the bubble surface along with the equilibrated bubble radius value. Fig. 2 shows the bubble growth and stability after the creation of the cavity. In the figure, only the O2 molecules are shown for the visualization purposes.

FIG. 2.

Bubble dynamics after the creation of the cavity in the system; A: t=200 ps; B: t=1 ns; C: t = 2 ns. The particular case shown corresponds to S = 200. The blue molecules represent O2.

FIG. 2.

Bubble dynamics after the creation of the cavity in the system; A: t=200 ps; B: t=1 ns; C: t = 2 ns. The particular case shown corresponds to S = 200. The blue molecules represent O2.

Close modal

The supersaturation value of the system changes with time because of the bubble growth. Since, the surface tension values were calculated based on the final bubble radius, only the final supersaturation values need to be considered. This is in contrast to the bubble experiment, where the amount of supersaturation remained constant because of the constant supply of gases. Four cases corresponding to the final supersaturation values of 1, 200, 400 and 850 were considered. Fig. 3 (a) shows that the relative surface tension of water decreases with an increase in the supersaturation ratio (final) or the concentration of the dissolved oxygen gas molecules. Svetovoy et al.6 observed that spontaneous combustion took place only for bubbles with a size smaller than 150 nm. For a bubble size of 150 nm, based on the calculated surface tension values from the present MD simulations, we predict the pressure of around 10 atm can be reached inside the nanobubbles (corresponding to 15 kHz and S = 200). For smaller bubbles, the internal pressure is even higher. Fig. 3 (b) shows how the pressure inside the bubble changes with the bubble diameter, for different supersaturation values. A big drop in the bubble pressure was obtained as the value of S was changed from 1 to 200 but after that the change in the bubble pressure with S is not significant.

FIG. 3.

a) Relative surface tension of water, γ(S)/γ(S=1, saturated water), as a function of the supersaturation; b) Internal bubble pressure as a function of the bubble diameter for different supersaturation values.

FIG. 3.

a) Relative surface tension of water, γ(S)/γ(S=1, saturated water), as a function of the supersaturation; b) Internal bubble pressure as a function of the bubble diameter for different supersaturation values.

Close modal

In summary, we performed molecular dynamics simulations to determine the surface tension of O2-supersaturated water. The surface tension of water was found to decrease with the amount of supersaturation, thus, the internal pressure in nanobubbles is much smaller than what would have been predicted using the planar-interface surface tension value of water. Nevertheless, the pressure inside nanobubbles is still very high, which may provide a suitable environment for spontaneous ignition and combustion to occur.

This work was supported by the National Science Foundation with Dr. Song-Charng Kong as the technical monitor.

1.
S.
Jain
,
O.
Yehia
, and
L.
Qiao
,
Journal of Applied Physics
119
,
094904
(
2016
).
2.
S.
Jain
,
W.
Park
,
Y. P.
Chen
, and
L.
Qiao
,
Journal of Applied Physics
120
,
174902
(
2016
).
3.
S.
Jain
,
G.
Mo
, and
L.
Qiao
,
Journal of Applied Physics
121
,
054902
(
2017
).
4.
S.
Lars
,
B.
Kevin
,
W.
Steffen
,
M.
Ewald
, and
R.
Paul
, in
39th Aerospace Sciences Meeting and Exhibit
(
AIAA
,
2001
).
5.
D. G.
Norton
and
D. G.
Vlachos
,
Chemical Engineering Science
58
,
4871
(
2003
).
6.
V. B.
Svetovoy
,
R. G. P.
Sanders
,
T. S. J.
Lammerink
, and
M. C.
Elwenspoek
,
Physical Review E
84
,
035302
(
2011
).
7.
V. B.
Svetovoy
,
R. G. P.
Sanders
, and
M. C.
Elwenspoek
,
Journal of Physics: Condensed Matter
25
,
184002
(
2013
).
8.
S.
Jain
,
A.
Mahmood
, and
L.
Qiao
, “
Quantifying heat produced during spontaneous combustion of H2/O2 nanobubbles
,”
2016 IEEE Sensor proceedings
,
Orlando, FL
,
Oct 30-Nov 3, 2016
.
9.
V.
Svetovoy
,
A.
Postnikov
,
I.
Uvarov
,
R.
Sanders
, and
G.
Krijnen
,
Energies
9
,
94
(
2016
).
10.
P.
Attard
,
The European Physical Journal Special Topics
,
1
(
2013
).
11.
B. D.
Johnson
and
R. C.
Cooke
,
Science
213
,
209
(
1981
).
12.
J.-Y.
Kim
,
M.-G.
Song
, and
J.-D.
Kim
,
Journal of Colloid and Interface Science
223
,
285
(
2000
).
13.
K.
Kikuchi
,
A.
Ioka
,
T.
Oku
,
Y.
Tanaka
,
Y.
Saihara
, and
Z.
Ogumi
,
Journal of Colloid and Interface Science
329
,
306
(
2009
).
14.
K.
Kikuchi
,
S.
Nagata
,
Y.
Tanaka
,
Y.
Saihara
, and
Z.
Ogumi
,
Journal of Electroanalytical Chemistry
600
,
303
(
2007
).
15.
K.
Kikuchi
,
H.
Takeda
,
B.
Rabolt
,
T.
Okaya
,
Z.
Ogumi
,
Y.
Saihara
, and
H.
Noguchi
,
Journal of Electroanalytical Chemistry
506
,
22
(
2001
).
16.
K.
Kikuchi
,
Y.
Tanaka
,
Y.
Saihara
,
M.
Maeda
,
M.
Kawamura
, and
Z.
Ogumi
,
Journal of Colloid and Interface Science
298
,
914
(
2006
).
17.
K.
Kikuchi
,
Y.
Tanaka
,
Y.
Saihara
, and
Z.
Ogumi
,
Electrochimica Acta
52
,
904
(
2006
).
18.
F.
Jin
,
J.
Li
,
X.
Ye
, and
C.
Wu
,
The Journal of Physical Chemistry B
111
,
11745
(
2007
).
19.
F.
Jin
,
J.
Ye
,
L.
Hong
,
H.
Lam
, and
C.
Wu
,
The Journal of Physical Chemistry B
111
,
2255
(
2007
).
20.
F.
Jin
,
X.
Ye
, and
C.
Wu
,
The Journal of Physical Chemistry B
111
,
13143
(
2007
).
21.
M.
Li
,
L.
Tonggu
,
X.
Zhan
,
T. L.
Mega
, and
L.
Wang
,
Langmuir
32
,
11111
(
2016
).
22.
H.
Kobayashi
,
S.
Maeda
,
M.
Kashiwa
, and
T.
Fujita
, in
Measurement and identification of ultrafine bubbles by resonant mass measurement method
,
2014
, p.
92320S
.
23.
J. L.
Parker
,
P. M.
Claesson
, and
P.
Attard
,
J. Phys. Chem.
98
,
8468
8480
(
1994
).
24.
D.
Andrienko
,
P.
Patricio
, and
O. I.
Vinogradova
,
J. Chem. Phys.
121
,
4414
4423
(
2004
).
25.
J.
Crassous
,
E.
Charlaix
,
H.
Gayvallet
, and
J. L.
Loubet
,
Langmuir
9
,
995
998
(
1993
).
26.
E.
Thormann
,
A. C.
Simonsen
,
P. L.
Hansen
, and
O. G.
Mouritsen
,
Langmuir
24
,
7278
7284
(
2008
).
27.
P.
Attard
,
Adv. Colloid Interface Sci.
104
,
75
91
(
2003
).
28.
M.
Holmberg
,
A.
Kühle
,
J.
Garnaes
,
K. A.
Mørc
, and
A.
Boisen
,
Langmuir
19
,
10510
(
2003
).
29.
A.
Simonsen
,
P.
Hansen
, and
B.
Klosgen
,
J. Colloid Interface Sci.
273
,
291
(
2004
).
30.
X. H.
Zhang
,
X. D.
Zhang
,
S. T.
Lou
,
Z. X.
Zhang
,
J. L.
Sun
, and
J.
Hu
,
Langmuir
20
,
3813
3815
(
2004
).
31.
X. H.
Zhang
,
G.
Li
,
N.
Maeda
, and
J.
Hu
,
Langmuir
22
,
9238
9243
(
2006
).
32.
X. H.
Zhang
,
N.
Maeda
, and
V. S. J.
Craig
,
Langmuir
22
,
5025
5035
(
2006
).
33.
S.
Yang
,
S. M.
Dammer
,
N.
Bremond
,
H. J. W.
Zandvliet
,
E. S.
Kooij
, and
D.
Lohse
,
Langmuir
23
,
7072
7077
(
2007
).
34.
N.
Ishida
,
M.
Sakamoto
,
M.
Miyahara
, and
K.
Higashitani
,
Langmuir
16
,
5681
5687
(
2000
).
35.
J. W. G.
Tyrrell
and
P.
Attard
,
Phys. Rev. Lett.
87
,
176104/1
176104/4
(
2001
).
36.
J. W. G.
Tyrrell
and
P.
Attard
,
Langmuir
18
,
160
167
(
2002
).
37.
M.
Alheshibri
,
J.
Qian
,
M.
Jehannin
, and
V. S. J.
Craig
,
Langmuir
32
,
11086
(
2016
).
38.
D.
Lohse
and
X.
Zhang
,
Reviews of Modern Physics
87
,
981
(
2015
).
39.
J. H.
Weijs
and
D.
Lohse
,
Physical Review Letters
110
,
054501
(
2013
).
40.
M. P.
Moody
and
P.
Attard
,
The Journal of Chemical Physics
120
,
1892
(
2004
).
41.
S.
He
and
P.
Attard
,
Physical Chemistry Chemical Physics
7
,
2928
(
2005
).
42.
S.
Plimpton
,
Journal of Computational Physics
117
,
1
(
1995
).
43.
H. J. C.
Berendsen
,
J. R.
Grigera
, and
T. P.
Straatsma
,
J. Phys. Chem.
91
(
24
),
6269
6271
(
1987
).
44.
J.-P.
Ryckaert
,
G.
Ciccotti
, and
H. J. C.
Berendsen
,
J of Comp Phys
23
,
327
341
(
1977
).
45.
H.
Andersen
,
J of Comp Phys
52
,
24
34
(
1983
).
46.
D.
Berthelot
,
Compt. Rendus
126
,
1703
1855
(
1898
).
47.
H. A.
Lorentz
,
Ann. Phys.
248
(
1
),
127
136
(
1881
).
48.
S.
Tanvir
,
S.
Jain
, and
L.
Qiao
,
Journal of Applied Physics
118
,
014902
(
2015
).
49.
J. H.
Irving
and
J. G.
Kirkwood
,
The Journal of Chemical Physics
18
,
817
(
1950
).
50.
W. G.
Hoover
,
Physical Review A
31
,
1695
(
1985
).
51.
S.
Nosé
,
The Journal of Chemical Physics
81
,
511
(
1984
).
52.
S.
Melchionna
,
G.
Ciccotti
, and
B.
Lee Holian
,
Molecular Physics
78
,
533
(
1993
).
53.
C.
Vega
and
E.
de Miguel
,
The Journal of Chemical Physics
126
,
154707
(
2007
).