Journal of Metallurgy

Journal of Metallurgy / 2011 / Article

Research Article | Open Access

Volume 2011 |Article ID 890321 |

Alexander L. Shimkevich, Inessa Yu. Shimkevich, "Molecular Dynamics Simulation of the Clustering of Minor Lead Additives in Liquid Sodium", Journal of Metallurgy, vol. 2011, Article ID 890321, 7 pages, 2011.

Molecular Dynamics Simulation of the Clustering of Minor Lead Additives in Liquid Sodium

Academic Editor: Seetharaman Sridhar
Received17 Dec 2010
Revised20 Mar 2011
Accepted15 Apr 2011
Published15 Jun 2011


A strong influence of minor lead additives on the liquid sodium microstructure is revealed in the molecular dynamics (MD) simulation of the Na0.98Pb0.02 alloy. The obtained results can be explained by the existence of lead-sodium clusters in liquid sodium built up by ionic bonds, Na+ , due to essential distinction of the alloy components in the electronegativity. On this reason, MD simulation of the Na0.98Pb0.02 alloy is carried out within the framework of a three-component bipolar model, Na + Na+ + , with Na Na+ recharging the nearest-neighbor particles of solvent in every 3 ps (an optimal period) during the numerical run.

1. Introduction

A nonideal solution of lead in liquid sodium, as an intermediate element in the alkaline-metal group, takes place in the whole composition range of Na-Pb alloy, but the strongest changes of thermodynamic properties are observed in the neighborhood of 20 and 50% at. of lead concentration [13]. It is connected with the presence of two kinds of clusters in the alloy: (Na4Pb)n and (NaPb)m, having an ionic bond, Na+–Pb, in diluted solutions and the covalent one, Na÷Pb, for the other compositions. They have explained this by the existence of Zintl’s cluster, (Na+)4(Pb4)4− [37]. The anionic tetrahedron, (Pb4)4−, is considered as a united particle with sodium cations, Na+, on the four faces of this tetrahedron.

However such anion is not in a steady state in liquid sodium [6, 7]. Therefore the application of Zintl’s cluster model [5] for understanding properties of the Na-Pb alloy is not enough. A new MD study as well as the experimental one of its microstructure and atomic dynamics by the method of neutron scattering is necessary. The knowledge about the Na-Pb alloy properties can be useful to develop a concept for improving the composition of the sodium coolant in particular by reducing its chemical activity in the environment and an automatic shut-down of sodium fires when sodium is used as a coolant of a fast nuclear reactor [8].

2. The Cluster Model of Na-Pb Alloy

It is known [9, 10] that the mixture of sodium and lead relates to strongly interacting systems with covalent bonds that are illustrated by a long homologous series of compounds, NanPbm, where n and m are integers. They influence the structural and thermodynamic properties of the binary alloy [17, 11, 12]. This effect can be displayed in analyzing the structural factors, SNN(q), SNC(q), and SCC(q), measured and calculated on various models [2, 5, 7]. It is also characterized by a sharp change of the alloy entropy [1] or an increase of the alloy electric resistance [3, 12] at ~ 0.25 as well as a deviation of alloy composition fluctuations from the ideal mixture shown in Figure 1 as a function of SCC(q)q→0 depending on the composition of the Na-Pb alloy [2].

One can see that the structural properties of Na-Pb alloy have distinctive features in the neighborhood of , which is equal to 0.2 and 0.5. This confirms the transitive character of sodium from easy lithium to heavy alkaline metals, Cs and K, and the existence of different clusters in the Na-Pb alloy such as (Na4Pb)k, (NaPb)l, and (NaPb4)m. In that case, the cluster model of solid spheres [2, 11] is the most realistic for considering the Na-Pb alloy in the whole range of composition. Moreover, the structural factor S(q) of Na0.8Pb0.2 alloy calculated on the solid sphere model of Takeda et al. [11] or Hafner et al. [9] agrees with the experimental data [11] obtained by scattering cold neutrons in samples of this alloy at 703 K (see Figure 2).

Thus, the structural factor of Na0.8Pb0.2 alloy has a subpeak in the neighborhood of q equal to 1.2 Å which is caused by spatial correlation of clusters in the alloy. The application of indirect Monte-Carlo method [14] confirms the presence of spatial microheterogeneity of Na0.8Pb0.2 alloy.

MD investigations of minor additives of lead in liquid sodium ( < 0.2) are a few. One can mention the calculated structural factor, S(q), of Na0.9Pb0.1 alloy obtained by Hafner et al. [9] and the experimental one obtained at neutron scattering [15] in the sodium alloy with 1.5 and 7.9% at. of lead.

At the same time, the character of heterogeneous structure of Na-Pb alloy remains unclear. Therefore we have applied a “partly ionized alloy” model [13, 16] for MD simulation of the Na-Pb alloy as lead anions and sodium cations in the liquid matrix of “neutral” sodium because of the approbation of this model proves to be successful for describing such systems.

3. MD Model of Na-Pb Alloy

There is a technique used for MD simulation of the Na1-x Pbx alloy of essentially different electronegative components which is based on the assumption that the additive one (Pb) is completely ionized at the electric neutrality of the system as a whole. This model describes a three-component system consisting of lead anions (Pb-), equal number of sodium cations (Na+), and neutral sodium atoms (Na) so that the total number N of particles in a MD cell is equal to + + .

Six pair potentials are necessary for describing the interaction of three types of particles: , , , , , and . We accept that is equal to which is the pair potential for the interaction of sodium particles in one-component system. For ionic pair potentials as , , and , the Born-Mayer-Huggins potential is used. It is applied earlier for simulating Na0.5Pb0.5 [5] and used here for describing the interaction of neutral sodium particle and lead anion; for a charge of sodium is equal to zero.

In the offered model, the charge of sodium cation can be “bound” rigidly to the particle (Na+) or handed to neighbor sodium atom by using the procedure of charge exchange among the nearest-neighbor particles ( ). This procedure can be a periodic or adaptive time function [17] for relaxing the charge perturbation of the system. This option allows approaching the three-component model of the Na-Pb alloy to a real electron exchange in the considered system.

4. Techniques for Numerical Runs

The MD simulation of Na1-x Pbx alloy (where x = 0.20, 0.09, 0.02) is carried out by MDMC code [17] in the frameworks of NVT ensemble of particles when all the lead atoms and the equal number of sodium atoms are oppositely ionized (Pb, Na+) so that the alloy remains neutral as a whole.

The periodic boundary conditions are used for a MD cube, and the interaction of particles is cut off in half of its edge. Ewald’s summation is used for calculating the ionic pair potentials in the approximation [5]. The equations of motion are solved by using Verlet’s algorithm in the velocity form [18] and a time step of 2·10−15 s. A numerical run is begun from random packing of particles in the MD cell, stopped at a given point of time for calculating required characteristics, and continued from the last atomic configuration. The temperature, pressure, kinetic and potential energy of the system, and their numerical error bars have been defined every 50 time steps. The root mean square deviations of energy and pressure did not exceed 1%.

The results are obtained by using Hasegawa’s potential [19] for = with a dielectric function [20], where RC is equal to 1.85 a.u. For the Born-Mayer-Huggins potentials: , , , and , the next ion diameters are used [5]: = 3.77 Å, = 3.53 Å, and = 3.70 Å, but = 3.20 Å is picked out by MD calculations. Obviously, the charge, , is equal to 1, = −1, and = 0.

The number of lead particles, , is defined by their atomic portion in the alloy as N, where N as a sum of + + is the total number of particles in the MD cell. The alloy of sodium and lead of 20 or 9% at. is simulated in the MD cell containing 5488 sodium particles ( + ). In the case of 2% at. of lead in liquid sodium, two systems are investigated—one with + = 9826 and the other with 18522 sodium particles.

Thus, we have had (1) N = 6860, = 4116, and = = 1372 in the MD cell for the Na0.8Pb0.2 alloy, (2) N = 6031, = 4945, and = = 543 for the Na0.91Pb0.09 one. Finally, we have investigated two MD ensembles for the Na0.98Pb0.02 alloy: (a) N = 10027, = 9625, and = = 201, (b) N = 18900, = 18144, and = = 378.

The Na0.8Pb0.2 alloy is simulated in MD cells with the molar volume of Ω ~ 247 a.u. [15] or Ω ~ 260 a.u., and the Na0.98Pb0.02 alloy is simulated in the MD cube with Ω ~ 297 a.u. (the edge length is equal to 76.13 Å).

The numerical experiment included three stages: (1)generating the random particle configuration at the beginning of the run, (2)putting the MD model of system in thermodynamic equilibrium, (3)calculating some characteristics of particle dynamics and fixing snapshots of atomic configurations of the system in points of time from 2 up to 100 ps in every 2 ps during the run.

Each run taken separately is carried out with or without the option of charge exchange (see above) every 3 or 4 ps. In the first case (with the option of charge exchange), the system is balanced every time as it is shown in Figure 3 which demonstrates a total energy change of the Na-Pb alloy after which the charge exchange is carried out.

It is visible that the MD simulation with the option of charge exchange of sodium particles in Na0.98Pb0.02 alloy deduces the system in the thermodynamic balance earlier and at higher level of its energy. It specifies the greater configuration disorder of system than it is without the option.

5. Results of MD Simulation of the Na-Pb Alloy

As marked above, there is the prepeak in the vicinity of q~ 1.2 Å in the structural factor, S(q), of the Na0.8Pb0.2 alloy (see Figure 2). This prepeak is caused by the spatial correlation of additive particles in liquid sodium. It proves also by molecular dynamics calculations in different models of the binary system (see Figure 4). It is visible that the results of molecular dynamics simulation of the Na-Pb alloy at 698 K in the model of partly ionized alloy [13] are closer to data of Hafner et al. [9].

On the contrary, such prepeak is absent in the S(q) curve of liquid sodium with the small concentration of lead, = 0.02, but at = 0.09, there is only a small bulge (see Figure 5). Its position in the model of partly ionized alloy [13] is in agreement with the same one of Hafner’s model [9], but a shift of the main S(q) peak to the left is more correspondent to the data on the neutron scattering [11] in sodium alloy with 1.5 and 7.9% at. of lead.

The 2% of lead atoms in sodium reduce diffusion mobility of sodium by 10% without changing the maximum position of vibration spectra. It specifies molecular bonds between ions in liquid sodium, the correlation of their Brownian motion, and a strong dissipation of oscillatory energy, for example, due to diffusive exchange of sodium particles in clusters (Na2Pb)n.

It is visible in Figure 3 that the thermodynamic relaxation kinetics of the particle charge exchange in every 3 and 4 ps is the same. Therefore, charge exchange periods less than 3 ps do not change the thermodynamic state of the microheterogeneous Na-Pb alloy.

As mentioned before, the MD simulation of the Na0.98Pb0.02 alloy with the option of recharging sodium particles is established earlier and at higher level of total energy (see Figure 3) and characterized by the greater configuration disorder of the system that is illustrated in Figures 610. These runs start with exactly the same initial state.

6. Discussion of the Obtained Results

One can see in Figures 7, 8, 9, and 10 that recharging sodium particles on simulating the Na0.98Pb0.02 alloy is the important option for adequate description of the clustering of the alloy, and the recharging period should not exceed 3 ps.

At the same time, the clustering of the ionic components (Na+, Pb) of the alloy illustrated in Figures 610 remains as sodium cations and lead anions quickly form chain structure of a percolation cluster which then gradually is transformed in separate cluster groups of different particle density. Here, it is important to note that recharging sodium particles (NaNa+) in every 3 ps and more often only reduces the cluster density but does not change the character of their reconstruction directed to the formation of steady-state colloidal particles (Na2Pb)n in a dilute solution of lead in liquid sodium.

These colloidal particles, as marked above, are characterized by molecular bonds between ions in liquid sodium and strong dissipation of oscillatory energy that is caused by an exchange of mobile sodium particles between clusters and the liquid matrix.

7. Conclusions

Using the three-component bipolar model, Na + Na++ Pb, with periodic recharging of the nearest-neighbor particles of solvent (NaNa+) has allowed us to carry out molecular-dynamic simulation of the Na0.98Pb0.02 alloy and to find (in the numerical runs) the strong influence of the minor additives of lead on the microstructure of liquid sodium.

The mechanism for clustering impurity particles in the dilute solution of lead in sodium is studied. It is shown that the essential distinction of electronegativity of the Na-Pb alloy components is the main reason for forming clusters on the base of ionic bonds, Na+–Pb.

Typical times of the colloid formation in any dilute solution of lead in sodium do not exceed 50 ps.


The authors are pleased to acknowledge Dr. A. S. Kolokol for discussing this work which is supported by the Russian Foundation of Basic Researches (Grant no. 09-08-00481a).


  1. K. Hoshino and W. H. Young, “On the entropy of mixing of the liquid Na-Pb alloy,” Journal of Physics F, vol. 11, no. 1, pp. L7–L9, 1981. View at: Publisher Site | Google Scholar
  2. S. Matsunaga, T. Ishiguro, and S. Tamaki, “Thermodynamic properties of liquid Na-Pb alloys,” Journal of Physics F, vol. 13, no. 3, pp. 587–595, 1983. View at: Publisher Site | Google Scholar
  3. W. van der Lugt, “Zintl ions as structural units in liquid alloys,” Physica Scripta, vol. 39, pp. 372–377, 1991. View at: Publisher Site | Google Scholar
  4. M. Tegze and J. Hafner, “Electronic structure of metallic and semiconducting alkali-metallead compounds,” Physical Review B, vol. 39, no. 12, pp. 8263–8274, 1989. View at: Publisher Site | Google Scholar
  5. H. T. J. Reijers, W. van Der Lugt, and M. L. Saboungi, “Molecular-dynamics study of liquid NaPb, KPb, RbPb, and CsPb alloys,” Physical Review B, vol. 42, no. 6, pp. 3395–3405, 1990. View at: Publisher Site | Google Scholar
  6. L. M. Molina, J. A. Alonso, and M. J. Stott, “Assembling alkali-lead solid compounds from clusters,” Journal of Chemical Physics, vol. 111, no. 15, pp. 7053–7061, 1999. View at: Google Scholar
  7. Y. Senda, F. Shimojo, and K. Hoshino, “The origin of the first sharp diffraction peak in liquid Na-Pb alloys: ab initio molecular-dynamics simulations,” Journal of Physics Condensed Matter, vol. 11, no. 10 A, pp. 2199–2210, 1999. View at: Google Scholar
  8. V. I. Subbotin, M. N. Arnol'dov, F. A. Kozlov, and A. L. Shimkevich, “Liquid-metal coolants for nuclear power,” Atomic Energy, vol. 92, no. 1, pp. 29–40, 2002. View at: Publisher Site | Google Scholar
  9. J. Hafner, A. Pasturel, and P. Hicter, “Simple model for the structure and thermodynamics of liquid alloys with strong chemical interactions. II. Chemical order and packing constraints,” Journal of Physics F, vol. 14, no. 10, pp. 2279–2295, 1984. View at: Publisher Site | Google Scholar
  10. Q. Wang, M. P. Iniguez, J. A. Alonso, and M. Silbert, “Charge transfer within clusters in liquid ionic alloys,” Journal of Physics: Condensed Matter, vol. 5, no. 26, pp. 4271–4282, 1993. View at: Publisher Site | Google Scholar
  11. S. Takeda, S. Harada, S. Tamaki, E. Matsubara, and Y. Waseda, “Structural study of liquid Na-Pb alloys by neutron diffraction,” Journal of the Physical Society of Japan, vol. 56, no. 11, pp. 3934–3940, 1987. View at: Google Scholar
  12. A. Thakur, N. S. Negi, and P. K. Ahluwalia, “Electrical resistivity of NaPb compound-forming liquid alloy using ab initio pseudopotentials,” Pramana—Journal of Physics, vol. 65, no. 2, pp. 349–358, 2005. View at: Google Scholar
  13. A. L. Shimkevich, I. Yu. Shimkevich, and V. V. Kuzin, “Molecular dynamics study of ternary system K+K+q+O-2q,” Condensed Matter Physics, vol. 2, no. 12, pp. 329–337, 1999. View at: Google Scholar
  14. S. Matsunaga, “Reverse monte carlo simulations in liquid Na-Pb and Na-Sn alloys,” Journal of the Physical Society of Japan, vol. 68, no. 7, pp. 2468–2469, 1999. View at: Google Scholar
  15. Frank Laboratory of Neutron Physics, Joint Institute for Nuclear Research. Annual report 2004,
  16. N. M. Blagoveshchenskii, V. A. Morozov, A. G. Novikov, V. V. Savostin, A. L. Shimkevich, and I. Y. Shimkevich, “Microscopic structure of liquid lead-potassium alloys: neutron-diffraction and molecular-dynamics investigation,” Physica B, vol. 364, no. 1–4, pp. 255–262, 2005. View at: Publisher Site | Google Scholar
  17. Yu. Shimkevich and A. L. Shimkevich, MDMC (Molecular-Dynamics and Monte-Carlo) software to study crystalline and irregular systems, Institute of Physics and Power Engineering, Preprint #2524, 1996.
  18. D. W. Heerman, Computer Simulation Methods in Theoretical Physics, Springer, 1986.
  19. M. Hasegawa, K. Hoshino, M. Watabe, and W. H. Young, “A new simple pseudopotential with applications to liquid metal structure factor calculations,” Journal of Non-Crystalline Solids, vol. 117-118, pp. 300–303, 1990. View at: Google Scholar
  20. S. Ichimaru and K. Utsumi, “Analytic expression for the dielectric screening function of strongly coupled electron liquids at metallic and lower densities,” Physical Review B, vol. 24, no. 12, pp. 7385–7388, 1981. View at: Publisher Site | Google Scholar

Copyright © 2011 Alexander L. Shimkevich and Inessa Yu. Shimkevich. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

More related articles

 PDF Download Citation Citation
 Download other formatsMore
 Order printed copiesOrder

Related articles

We are experiencing issues with article search and journal table of contents. We are working on a fix as to remediate it and apologise for the inconvenience.

Article of the Year Award: Outstanding research contributions of 2020, as selected by our Chief Editors. Read the winning articles.