1. Introduction: phenomenology of laser ablation in liquid
Laser ablation of solids in liquids is considered a very efficient technique for the synthesis of nanocrystals. The most important feature is the extreme confinement effect of liquid, i.e. the liquid restricts the expansion of the plasma plume. The extreme confined conditions and induced high-pressure region favour the formation of unusual metastable phases (Wang & Yang Reference Wang, Yang and Yang2012). These advantages allow the designer to combine selected solid targets and liquid to fabricate compound nanostructures with desired functions.
The different phases involved in the process can be summarized as follows (see figure 1). At the initial stage of the interaction of the high-energy laser with the interface between the solid and the liquid, a plasma plume is generated and due to its high kinetic energy, it expands. In liquid-phase laser ablation, the liquid environment presents a strong confinement effect limiting the expansion rate. The result of this effect is a very high number density inside the plasma (with
$n_{pl}>10^{25}~\text{m}^{3}$
), coupled with a fast decrease of the temperature due to the predominance of recombination phenomena, and a persistence of plasma in the plume for a few hundred nanoseconds, which is a lifetime 10 times shorter than that of plasma in air (Wang & Yang Reference Wang, Yang and Yang2012). Another subsequent effect of the localized high thermal exchange between plasma and liquid is that a small amount of the surrounding liquid is vaporized to form a cavitation bubble between the plasma front and the liquid at about
$1~{\rm\mu}\text{s}$
. The bubble expands to the maximum size, then it shrinks and collapses after about
$300~{\rm\mu}\text{s}$
.
Sketch of the main stages involved in a laser-induced plasma by target ablation in liquid.
Although the laser-induced cavitation bubble plays an important role in the nanoparticle dynamics (during the bubble shrinking phase nanoparticles can acquire compactness), it is now well established (Sasaki & Takada Reference Sasaki and Takada2010; Soliman, Takada & Sasaki Reference Soliman, Takada and Sasaki2010) that the formation of nanoparticles takes place long before, during the plasma-phase expansion (in addition it is possible that some clusters have issued directly from solid material but this formation mechanism is more typical of pico- and femto-second lasers). Therefore one can distinguish three different phases in the nanoparticle synthesis (Amendola & Meneghetti Reference Amendola and Meneghetti2013):
-
(a) Nucleation condensation around ion-seeds in the initial strongly coupled plasma ablated from the target and confined by the liquid environment. This process leads to the formation of liquid-phase embryos of few atoms (the system is known as ‘misty plasma’ (Coppins Reference Coppins2010)). This process is analysed in § 2.
-
(b) Growth of the drop mediated by the surrounding plasma medium.
-
(c) Crystallization of the drop passing through the cavitation bubble boundary.
The nucleation and growth of the nanoparticle can be described through a thermodynamic approach (classical theory of nucleation) or by a kinetic approach (Debenedetti Reference Debenedetti1996). The first assumes that the embryos are in equilibrium to calculate the condensation and evaporation coefficients and that they behave as bulk material (approximation of capillarity), while the second approach is more rigorous and is based on first principles.
In this work, the second approach has been considered, developing a Monte Carlo method for the condensation and evaporation of nanoparticles in plasma taking into account the electrostatic attraction–repulsion driving mechanism due to the charging of nanoparticles. In § 2 the nucleation condensation theory will be recalled stressing the role of the plasma medium, while in § 3 the kinetic of growth will be presented in detail, showing the condensation mechanism in plasma and the evaporation rate as resulting from Reiss–Kegel–Katz (RKK) theory (Reiss, Kegel & Katz Reference Reiss, Kegel and Katz1997). The numerical Monte Carlo algorithm used is described in § 4, while finally in § 5 results of the numerical models applied to the case of a nickel target ablated in liquid using a nanosecond laser pulse will be presented and discussed.
2. Particle nucleation condensation in plasma
2.1. Homogeneous nucleation
The theory is based on the assumptions that the new-phase embryos are perfectly spherical, non-deformable and not interacting with each other. The phase transition is driven by the work (difference between the cohesion energy of atoms in the new condensed phase and the energy barrier due to surface tension) required to produce a nucleus of size
$r$
composed of
$N$
atoms:
where
$A=4{\rm\pi}r^{2}$
is the surface area of the embryo,
${\it\gamma}$
its surface tension and
${\rm\Delta}{\it\mu}$
the chemical potential variation between the gas and the liquid phase of the embryo (
${\rm\Delta}{\it\mu}=kT\ln S$
where
$S=P_{v}/P_{ev}$
is the supersaturation ratio, considering the liquid phase to be at the equilibrium vapour pressure). If each atom of the liquid embryo occupies a volume
$V$
(
$V=1/{\it\varrho}$
, where
${\it\varrho}$
is the density of the embryo), we can express the work as a function of the radius:
The critical size of nucleation, the critical energy, the transition probability and the rate of nucleation are all determined by (2.2). The critical size
$r^{\ast }$
is calculated by imposing
$\partial {\rm\Delta}G(r)/\partial r=0$
leading to:
$$\begin{eqnarray}r^{\ast }=\frac{2{\it\gamma}V}{kT\ln \left(\displaystyle \frac{P_{v}}{P_{ev}}\right)}.\end{eqnarray}$$
The critical radius implies that larger embryos grow spontaneously, while below this value they remain spontaneous density fluctuations without condensing in the new phase, but with a natural tendency to disintegrate. Note that the critical radius is inversely proportional to
${\rm\Delta}{\it\mu}$
, i.e. the more metastable the plasma is, the smaller the critical nuclei are. The size distribution of embryos formed follows the Maxwell–Boltzmann behaviour:
2.2. Ion-induced nucleation
A mechanism for the role of ions in condensation has been proposed by previous investigators (Russell Reference Russell1969) and it has been associated to the formation of clouds in the upper atmosphere (Keesee Reference Keesee1989). The presence of plasma as a medium favours the nucleation because it lowers the work required (2.2) by adding a term taking into account the presence of ionization:
where
${\it\varepsilon}_{c}$
is the dielectric constant of the clusters (nickel in the present case),
$Q$
is the ion charge,
$r_{a}$
is the radius of a jacketed ion and
${\it\varepsilon}=1$
has been assumed for the vapour phase.
As a consequence the critical radius
$r^{\ast }$
(Girshick Reference Girshick1994; Tillack, Blair & Harilal Reference Tillack, Blair and Harilal2004) is reduced and the nucleation rate increases. At sufficiently high saturation ratios, the ion-induced nucleation can becomes ‘barrierless’, and every ion is potentially a seed site.
3. Kinetic theory of particle growth
The kinetic growth framework is based on attachment and detachment of single atoms or of clusters of atoms, according to:
Here,
$N_{j}$
denotes a particle of size
$j$
, with
$p_{\mathit{agg}}(j,k)$
the rate constant of the two-particle attachment reaction of a particle of size
$j$
to a particle of size
$k$
(reaction from left to right), and
$p_{ev}(j+k,j)$
the rate constant of the one-particle detachment of a particle of size
$j$
from a particle of size
$j+k$
(reaction from right to left). The reaction rates often used are based on short-range collisional interaction, while in our case, an electrostatic attraction between the negative-charged nanoparticle and the plasma ions is the driving mechanism.
The kinetic approach is based on the embryo growth through the explicit sequential adhesion of ions, knowing the probability of adhesion and evaporation as the nanoparticle grows and as the plasma parameters (density and temperature) in which it is immersed change. The nanoparticle growth occurs through the coupling between the granule charge and the dynamics of the plasma plume emitted from the target and confined by the liquid medium. In fact, the charging time of a 1 nm size embryo
${\it\tau}_{ch}$
(Cui & Goree Reference Cui and Goree1994) is comparable with the plasma recombination (the dominant plasma recombination channel is represented by a three-body recombination
$2e+Ni^{+}\rightarrow e+Ni$
with
${\it\tau}_{\mathit{rec}}=2\times 10^{-12}~\text{s}$
for a plasma with density
$n_{pl}=10^{18}~\text{cm}^{-3}$
and temperature
$T_{pl}=6000~\text{K}$
). In the present model, the plasma density is assumed to exponentially decay with time:
starting from an initial density
$n_{pl,\mathit{max}}$
, while the temperature can undergo an isothermal (
${\it\gamma}=1$
) or an adiabatic (
${\it\gamma}=5/3$
) cooling expansion:
where
$T_{pl}(0)$
is the initial temperature. The temperature can undergo an isothermal expansion in the case of a very fast heat exchange between the plasma plume front and the liquid environment. In the opposite case, for a very slow energy exchange between the plasma particles and the liquid (for example in the case of a liquid characterized by a low thermal conductivity), the expansion is adiabatic.
A strong assumption is made in this first model: plasma is considered homogeneous and both electron and ion components are considered in equilibrium at temperature
$T_{pl}$
. A more sophisticated model may couple the nanoparticle growth with the simulation of the expansion of the plasma plume confined by the liquid. Such models are completely missing up to now and the only attempts to simulate the expansion of the plasma plume have been made through the gas medium at various pressures (Itina et al.
Reference Itina, Hermann, Delaporte and Sentis2002; Itina, Sentis & Marine Reference Itina, Sentis and Marine2006; Harilal et al.
Reference Harilal, Miloshevsky, Diwakar, LaHaye and Hassanein2012).
The negative charge of the nucleus (due to the much larger electron mobility) attracts ions around it, allowing an agglomerative (single positive ion) type of growth in place of an aggregative (negatively charged cluster) type (see figure 2):
-
(i)
$N_{j}+N_{1}\rightarrow N_{j+1}$
(agglomeration).
The charge on a particle is closely related to its size and greatly depends on the local plasma conditions (Kundrapu et al.
Reference Kundrapu, Li, Shashurin and Keidar2012). The confined plasma ensures a high degree of plasma uniformity which is reflected in a good quality nanoparticle size distribution function. Every particle attains a certain floating potential that ensures equality of the electron and positive ion currents towards the nanoparticle surface. For a spherical nanoparticle with radius
$r_{p}$
, these currents can be calculated by means of the orbital-motion-limited (OML) theory when the condition
$r_{p}\ll {\it\lambda}_{D}$
applies (
${\it\lambda}_{D}$
represents the Debye length). For the collection of Maxwellian electrons and ions by a negatively charged particle, the orbit-limited currents are given by (Cui & Goree Reference Cui and Goree1994; Lampe Reference Lampe2001)
where
$C$
is the capacitance of a spherical object in plasma.
Sketch of electron and ion plasma particle collection on the nanoparticle surface.
Finally, the rate of evaporation from the nucleus has to be considered:
-
(ii)
$N_{j+1}\rightarrow N_{j}$
(evaporation)
through the classical RRK theory (Itina et al.
Reference Itina, Sentis and Marine2006, Reference Itina, Gouriet, Zhigilei, Noel, Hermann and Sentis2007; Gouriet, Sentis & Itina Reference Gouriet, Sentis and Itina2009; Jansen et al.
Reference Jansen, Wysong, Gimelshein, Zeifman and Buck2010). The evaporation rate of a cluster with
$N$
atoms is given as
where
$p_{0}$
is the typical molecular vibrational frequency,
$s_{N}$
is the number of vibrational degrees of freedom (1 for dimers and
$3N-6$
for larger clusters),
$g_{N}$
is the number of surface molecules,
$E_{N}$
is the absolute binding energy of the cluster and
$D_{N}$
is the molecular evaporation energy. In the present case, molecular dynamics results (Gouriet et al.
Reference Gouriet, Sentis and Itina2009) have been used with the following fitting formula as a function of the number of atoms
$N$
in the nickel cluster, and as a function of the plasma parameters, density and temperature:
where
$T_{pl}$
and
$n_{pl}$
are in units of K and
$\text{m}^{-3}$
, respectively.
4. Monte Carlo numerical model
The Monte Carlo model developed is based on the OML theory described in the previous section and it consists of an evolution of numerical models used to study electrostatic probes and dust in plasma (Taccogna, Longo & Capitelli Reference Taccogna, Longo and Capitelli2003, Reference Taccogna, Longo and Capitelli2004; Taccogna Reference Taccogna2012; Taccogna & Mizzi Reference Taccogna and Mizzi2014). The algorithm involves the self-consistent change of charge, surface electric potential and radius of the granule during the expansion of the plasma. For every simulated nanoparticle
$p$
(starting as an embryo of
$N_{p}$
atoms with charge
$Q_{p}=0$
and surface potential
${\it\phi}_{p}=0$
), the Monte Carlo algorithm consists of the following instructions during the time step
${\rm\Delta}t$
(see figure 3):
-
(a) calculate the total plasma particle (electron/ion) collisional rate (probability per unit time) with the granule:
$p_{\mathit{agg}}=p_{e}+p_{i}=I_{e}/e+I_{i}/e$
using (3.4)–(3.5), and of an atom evaporation rate from the granule
$p_{ev}$
using (3.8); -
(b) if a plasma particle attachment happens, choose the collisional partner, electron or ion:
$p_{{\it\alpha}}/p_{\mathit{agg}}$
(
${\it\alpha}=e,i$
), and make the following changes to the nanoparticle:-
(i) the nanoparticle charge:
$Q_{p}=Q_{p}\mp e$
(
$-$
for electron and
$+$
for ion; in the latter case update the number of ions in the granule
$N_{p}=N_{p}+1$
); -
(ii) the radius of the nucleus, which changes with the number of ions in the nanosphere:
$r_{p}=(aN_{p})^{1/3}$
, where
$a=(3\times 10^{-3}A)/(4{\rm\pi}{\it\rho}N_{A})$
with
$A$
the atomic number,
${\it\rho}$
the density (
$\text{g cm}^{-3}$
) and
$N_{A}$
Avogadro’s number; -
(iii) the surface electric potential as prescribed by (3.6):
${\it\phi}_{p}=Q_{p}/(4{\rm\pi}{\it\varepsilon}_{0}r_{p}(1+r_{p}/{\it\lambda}_{D}))$
;
-
-
(c) in the case of evaporation update the number of atoms in the cluster
$N_{p}=N_{p}-1$
and thus the radius and the surface electrostatic potential as at point (b); -
(d) update the plasma density and temperature following (3.2) and (3.3), respectively.
Flow chart of the Monte Carlo numerical model for nanoparticle growth in plasma.
The simulation for every nanoparticle advances using a time step that is not fixed. It is chosen so that it corresponds to the time interval between collection of plasma particles or evaporation of an atom from the cluster consistent with the probability
$p_{\mathit{tot}}=p_{\mathit{agg}}+p_{ev}$
:
where
$\mathit{Rand}$
is a random number uniformly distributed between 0 and 1. Once the plasma is extinguished, the plasma-induced agglomeration and evaporation probability goes to zero (the iteration of the algorithm stops) and the nanoparticles are driven through thermophoresis forces towards the plasma–liquid interface where a steep temperature drop makes them crystallize and be released into the newly forming bubble (Girshick Reference Girshick1994).
5. Results and discussion
The parameters used for the simulation are those typical of a nanosecond laser pulse on a nickel target in water (irradiance
$I=140~\text{GW cm}^{-2}$
, laser pulse duration
${\it\tau}=6~\text{ns}$
) (De Giacomo et al.
Reference De Giacomo, Dell’Aglio, Santagat, Gaudiuso, De Pascale, Wagener, Messina, Compagnini and Barcikowski2013; Dell’Aglio et al.
Reference Dell’Aglio, Gaudiuso, De Pascale and De Giacomo2015). The plasma density reaches its maximum (
$n_{pl,\mathit{max}}=2\times 10^{25}~\text{m}^{-3}$
) after
$t_{\mathit{max}}=20~\text{ns}$
, while it is extinguished in
${\it\tau}=100~\text{ns}$
. The initial temperature is set at
$T_{pl,0}=6000~\text{K}$
.
Figure 4 shows the behaviour of plasma density (black line) and temperature (red lines) following (3.2) and (3.3) using the above parameters in the case of isothermal (full red line) and adiabatic (broken red line) conditions.
Temporal evolution of plasma density (black line) and temperature (red lines) in the plume using (3.2) and (3.3) with
$n_{pl,\mathit{max}}=2\times 10^{25}~\text{m}^{-3}$
),
$t_{\mathit{max}}=20~\text{ns}$
,
${\it\tau}=100~\text{ns}$
,
$T_{pl,0}=6000~\text{K}$
and
${\it\gamma}=1$
(isothermal, full line) and
${\it\gamma}=5/3$
(broken line).
The kinetic growth model shows a weak dependence on the starting size of embryos. Whatever the number of atoms in the initial embryo (between 1 and 500 atoms, i.e. for an initial radius smaller than 1 nm) all the nanoparticles reach the same final radius.
The typical growth is shown in figure 5 where the temporal evolution of the nanoparticle radius is reported starting with an initial radius of
$r_{p,0}=0.6~\text{nm}$
(
$N=100$
atoms). The black line refers to the case of isothermal plasma expansion while the red line has been obtained using an adiabatic cooling expansion. It is evident that for both the cases the growth starts at about 0.5 ns after the appearance of the embryo and it stops at 300 ns. At this point the plasma no longer has a role in driving the agglomeration and can be considered extinguished regarding nanoparticle growth. While in the isothermal case the nanoparticle radius reaches a value of about
$r_{p}=4~\text{nm}$
, for the adiabatic expansion the nanoparticle is bigger having a radius
$r_{p}=7~\text{nm}$
(in De Giacomo et al.
Reference De Giacomo, Dell’Aglio, Santagat, Gaudiuso, De Pascale, Wagener, Messina, Compagnini and Barcikowski2013; Dell’Aglio et al.
Reference Dell’Aglio, Gaudiuso, De Pascale and De Giacomo2015 nanoparticle radii in the range of 5–10 nm have been measured). The two curves diverge at about 40 ns, when the adiabatic case is characterized by a lower plasma temperature than the isothermal case. The lower temperature inhibits the evaporation that controls the nanoparticle size. In such a way condensation remains the dominant process allowing the synthesis of larger nanospheres.
Temporal evolution of the radius of a nickel nanoparticle starting with an initial radius of
$r_{p,0}=0.6~\text{nm}$
for an isothermal (black line) and an adiabatic (red line) plasma expansion.

In order to analyse further in detail the growth mechanism, in figure 6 the ion agglomeration and atom evaporation rates are reported as a function of time for the isothermal case. For
$t\leqslant 1~\text{ns}$
both the rates are irrelevant, so preventing any electrostatic-induced growth. Later, and up to 200 ns, the ion collection rate is larger than the atom evaporation allowing a continuous growth until the agglomeration rate, having reached its maximum, starts to drop and becomes equal to the evaporation rate. At this point the plasma quenching stops any further agglomeration and evaporation from the surface and the nanoparticle radius reaches a steady value.
Temporal evolution of the ion agglomeration (black line) and atom evaporation (red line) rate from the nickel nanoparticle whose evolution is reported in figure 5 (isothermal case).

As pointed out before in § 3, the electrostatic attraction between the negatively charged nanoparticle and the positive plasma ions is the driving growth mechanism and this is evident from figure 7 where both the nanoparticle charge (black curve) and its surface potential (red curve) are reported as a function of time for the isothermal plume expansion case. In the first 100 ns, the nucleated embryo is negatively charged due to the electron collection on its surface. As a consequence the surface potential becomes negative,
${\it\phi}_{p}=-2.5~\text{V}$
compared to the plasma reference potential. From this instant, the nanoparticle starts to collect ions and aggregations takes place, leading to a small increase of the nanoparticle negative charge, that nevertheless maintains a steady-state negative value of
$Q_{p}=-8e$
that gives the nanoparticle a high stability.
Temporal evolution of the nanoparticle charge
$Q_{p}$
(black line) and surface electric potential
${\it\phi}_{p}$
(red line) of the nickel nanoparticle whose evolution is reported in figure 5 (isothermal case).
6. Conclusions
In this work, a kinetic model of the nanoparticle growth in a laser-induced plasma by target ablation in liquid has been developed by means of a Monte Carlo algorithm. The mechanism involved in the agglomeration consists of the nanoparticle charging by electron collection, with subsequent ion attraction to maintain the zero current collection on the nanosphere surface. In addition, an evaporation mechanism is taken into account by means of a molecular dynamics fitting result. Both the competing processes are a function of the plasma parameters, density and temperature, whose evolution fits experimental measurements. The model highlights the temporal evolution of the nanoparticle size, showing the important effect of the confining role of the liquid on the plasma plume expansion (more isothermal than adiabatic) in creating small nanoparticles with an equilibrium negative charge and surface potential. Future works will be devoted to the inclusion of a plasma plume model in order to increase the level of self-consistency of the kinetic nanoparticle growth.
Acknowledgements
This research has received funding from the Italian Ministry of Education, University and Research (MIUR) under the project PRIN 2010ERFK 007. The authors would like to acknowledge the LIP group of CNR-NanoTec (A. De Giacomo, M. Dell’Aglio, R. Gaudiuso and G. Valenza) for fruitful discussions.

