Thermodynamics of Point Defects in Solids and Relation with the Bulk Properties: Recent Results

Abstract

For several decades, the crucial question has arisen as to whether there exists any direct interconnection between the thermodynamic parameters of point defects in solids with the bulk properties of the solid under investigation. To answer this important question, an interrelation of the defect Gibbs energy $g^i$ in solids with bulk properties has been proposed almost half a century ago. Considering that $g^i$ corresponds to an isobaric and isothermal process, this interrelation states that, for different processes (defect formation, self-diffusion activation, and heterodiffusion), $g^i$ is proportional to the isothermal bulk modulus B and the mean volume per atom $\mathrm{\Omega }$, termed cB$\mathrm{\Omega }$ model. Here, we review several challenging applications of this interrelation that appeared during the last decade (2011–2021), including high pressure diamond anvil measurements, high T${}_c$ superconductors, nuclear fuels, and materials for micro-electronics devices, applications of usefulness in Geophysics and Seismology, a problem of major technological interest, search for compositions of better target properties in Cu-Co-Si alloys via machine learning as well as two independent studies on the physical origin of this interrelation that has been further strengthened during the last few years.

1. Introduction

The question of whether the defect Gibbs energy $g^i$ (where i denotes the corresponding process, i.e., defect formation, f, migration, m, or self-diffusion activation, $act$) is explicitly interrelated with bulk properties in solids is of paramount importance in the studies of crystals. Almost half a century ago, a model [1,2,3,4,5,6,7,8,9,10,11]—hereafter called cB$\mathrm{\Omega }$ model—argued (e.g., see Ref. [12]) that $g^i$ should be proportional to the isothermal bulk modulus B and the mean volume $\mathrm{\Omega }$ per atom, which differs essentially [12] (for a review, see also Ref. [13]) from an earlier model by Zener [14,15,16] in which $g^i$ is related with the shear modulus, $\mu$.

Following a strict foundation of thermodynamics of point defects in solids [12], the defect formation parameters are defined by comparing a real (i.e., containing defects) crystal either to an isobaric ideal (i.e., not containing defects) crystal or to an isochoric ideal crystal. Thus, there are two different families of defect formation parameters, which are interconnected through thermodynamic relations in which the isothermal bulk modulus B plays an important role [12]. The isobaric parameters are defined in terms of the corresponding Gibbs energy ($g^f$) as follows:

$$s^f=-{\displaystyle \frac{d{g}^f}{dT}}{|}_P$$

(1)

$$h^f=g^f-T{\displaystyle \frac{d{g}^f}{dT}}{|}_P$$

(2)

and hence $h^f=g^f+T{s}^f$,

$${\upsilon }^f={\displaystyle \frac{d{g}^f}{dP}}{|}_T$$

(3)

where (P and T denote the pressure and temperature, respectively) and $s^f$, $h^f$ and ${\upsilon }^f$ designate the defect formation entropy, enthalpy, and volume, respectively.

According to the cB$\mathrm{\Omega }$ model, the defect Gibbs energy $g^i$ is interconnected with the bulk properties of the solid through the relation [2,3,4,5,6,7,8,9,10,11] (cf., for its justification, see Ref. [12]):

$$g^i=c^{i}B\mathrm{\Omega }$$

(4)

where $c^i$ is dimensionless, which—to the first approximation—can be considered as independent of temperature and pressure (cf. this approximation holds if the ratio $\Delta P/B$—where $\Delta P$ refers to the pressure range under consideration—is appreciably smaller than unity; otherwise, a correction factor should be also taken into account, see p. 126 of Ref. [12]). The superscript i in Equation (4) refers, as mentioned, to the defect processes under consideration, i.e., f, $act$, m (formation, activation and migration, respectively), and of course $c^f\ne c^m\ne c^{act}$. By inserting Equation (4) into Equations (1)–(3), we find

$$s^i=-c^i\mathrm{\Omega }(\beta B+{\displaystyle \frac{dB}{dT}}{|}_P)$$

(5)

$$h^i=c^i\mathrm{\Omega }(B-T\beta B-T{\displaystyle \frac{dB}{dT}}{|}_P)$$

(6)

$${\upsilon }^i=c^i\mathrm{\Omega }({\displaystyle \frac{dB}{dP}}{|}_T-1)$$

(7)

where $\beta$ is the thermal (volume) expansion coefficient of the bulk material. These relations suggest that the quantities $s^i/h^i$, ${\upsilon }^i/h^i$ and ${\upsilon }^i/g^i$ solely depend on bulk properties. Furthermore, we note that, for the defect volume ${\upsilon }^i$, its thermal expansion coefficient ${\beta }^i=\frac{1}{{\upsilon }^i}{\displaystyle \frac{\partial {\upsilon }^i}{\partial T}}{|}_P$ and compressibility ${\kappa }^i=-\frac{1}{{\upsilon }^i}{\displaystyle \frac{\partial {\upsilon }^i}{\partial P}}{|}_T$ can be explicitly calculated from Equation (7) and—as shown in detail in Ref. [12]—exceed markedly the corresponding quantities $\beta$ and $\kappa$, respectively of the bulk solid. They are given by

$${\kappa }^{act}=\kappa -\frac{\frac{d^{2}B}{d{P}^2}{|}_T}{\frac{dB}{dP}{|}_T-1}$$

(8)

and

$${\beta }^{act}=\beta +\frac{\frac{d}{dT}\left\{\frac{dB}{dP}{|}_T\right\}}{\frac{dB}{dP}{|}_T-1},$$

(9)

which can be alternatively written as:

$$\frac{{\kappa }^i}{\kappa }=1-\frac{B\frac{d^{2}B}{d{P}^2}{|}_T}{\frac{dB}{dP}{|}_T-1}$$

(10)

and

$$\frac{{\beta }^i}{\beta }=1+\frac{\frac{d}{dT}\left\{\frac{dB}{dlP}{|}_T\right\}}{\beta \frac{dB}{dP}{|}_T-1}.$$

(11)

When a single mechanism is operating in mono-atomic crystals, the self-diffusion activation process can be expressed through the activation Gibbs energy $g^{act}$, which stands for the sum of the Gibbs energies $g^f$ and $g^m$ corresponding to the formation and migration processes, respectively. The activation entropy $s^{act}$ and the activation enthalpy $h^{act}$ are then defined accordingly (e.g., see Ref. [12]): $s^{act}=-{\displaystyle \frac{d{g}^{act}}{dT}}{|}_P$ and $h^{act}=g^{act}+T{s}^{act}$ and the diffusion coefficient D is given by

$$D=\gamma {\alpha }^2\nu exp(-{\displaystyle \frac{g^{act}}{k_{B}T}})$$

(12)

and, inserting Equation (4) into Equation (12), we get

$$D=\gamma {\alpha }^2\nu exp(-{\displaystyle \frac{c^{act}B\mathrm{\Omega }}{k_{B}T}})$$

(13)

where $\gamma$ is a numerical constant depending on the diffusion mechanism and the structure, $\alpha$ denotes the lattice constant, and $\nu$ the attempt frequency which is of the order of the Debye frequency ${\nu }_D$. The plots $lnD$ versus $1/T$ (for P = const) and $lnD$ versus P (for T = const) are usually found to be straight lines that reflect that $g^{act}$ varies linearly upon increasing the temperature (for P = const) and pressure (for T = const), respectively, but in some cases, e.g., Na, these plots deviate strongly from linearity.

It is the scope of this feature article to present the most striking applications of the cB$\mathrm{\Omega }$ model that appeared in the literature in various scientific disciplines during the last decade, as well as two studies on the foundation of the physical origin of the cB$\mathrm{\Omega }$ model which have just been published during the last two years; see the block diagram in Figure 1. In view of the multitude of the scientific disciplines discussed in each box of the block diagram, we insert for reader’s convenience the section at which each of the matters is described.

2. Materials and Methods

The temperature variation of D at a constant pressure can be deduced in the frame of the cB$\mathrm{\Omega }$ model as follows (for a review, see also Ref. [13]): When the value $D_1$ is determined experimentally at a temperature $T_1$, the value of $c^{act}$ can be estimated because the factor $\gamma {\alpha }^2\nu$ is approximately known. An error of around 2 when setting $\nu$ equal to ${\nu }_D$ does not change the value of $c^{act}$ significantly [8]. Hence, once the value of $c^{act}$ has been determined from $D_1$, the value of $D_2$ for any temperature $T_2$ can be estimated from Equation (13) if the elastic and expansivity data are known at this temperature. Once $c^{act}$ has been estimated $D_1$, the values of $s^{act}$, $h^{act}$, and ${\upsilon }^{act}$ can be computed at any temperature from Equations (5)–(7), respectively. The D value can be studied at any pressure for $T=$ const when $c^{act}$ is estimated through Equation (13) from $D_1$ measured at a single pressure $P_1$.

3. Applying cB$\mathrm{\Omega }$ Model to High Pressure Diamond Anvil Measurements

Lithium fluoride (LiF) is a unique crystal that exhibits the largest band gap. It is believed that it remains transparent to visible light under stresses in excess of 1000 GPa [17]. Before 2014, the high-pressure static compression behavior of LiF has been determined only up to 9 GPa under hydrostatic conditions [18]. In 2014, the compression behavior of LiF and its equation of state have been determined to 92 GPa by X-ray diffraction in a diamond anvil cell [19]. Recently [20], an independent study reported the equation of state of several solids (including LiF, NaCl, and Pb) measured under ambient temperature in the Mbar range using diamond anvil cells with a quasi-hydrostatic pressure transmitting medium (helium or neon) and coupled with X-ray diffraction.

In 2020, the applicability of the cB$\mathrm{\Omega }$ model was studied by using the equation of state of LiF that has been obtained, as mentioned above, in Ref. [20]. Equation (7) can be rewritten as:

$$v^i=\frac{g^i}{B}\left.{\left(\frac{dB}{dP}\right|}_T-1\right)$$

(14)

or

$$\frac{v^i}{g^i}=\frac{1}{B}\left.{\left(\frac{dB}{dP}\right|}_T-1\right)$$

(15)

The Rydberg–Vinet equation of state of LiF in the Mbar range that was reported in the recent work by Dewaele [20] (third column of his Table 1) gives at T = 300 K: B = 62.3 GPa and $\frac{dB}{dP}{|}_T=5.01$ ($\mathrm{\Omega }$ = 16.391 Å${}^3$). In a previous work, Dong et al. [19], after compression of LiF to 92 GPa, found B = 66.2 GPa and $\frac{dB}{dP}{|}_T=4.6\pm 0.1$ through the 3rd-order Birch–Murnaghan equation of state.

McKeever and Lilley [21] measured the thermally stimulated polarisation currents (TSPC) and the thermally stimulated depolarisation currents (TSDC) of LiF doped with Mg${}^{2+}$. Their TSPC and TSDC curves exhibited the maximum (approximately around at T${}_{max}$ = 211 ± 3 K) corresponding to a single relaxation process, the analysis of which resulted in [21]: $h_{bv}^m$ = 0.66 eV and ${\tau }_0^{-1}=1.9\times {10}^{14} {\mathrm{s}}^{-1}$ of the parameters ${\tau }_0$, $h_{bv}^m$ in the Arrhenius relation:

$$\tau ={\tau }_0\exp (h_{bv}^m/kT)$$

(16)

where ${\tau }_0$ denotes the pre-exponential factor and $h_{bv}^m$ denotes the enthalpy for the bound cation vacancy motion.

Following Ref. [22], the quantity $s_{bv}^m$ for the bound vacancy motion in the (re)orientation process of the dipole ‘Mg${}^{2+}$-cation vacancy’ is estimated from

$$s_{bv}^m=k\ln ({\tau }_0^{-1}/2\nu )$$

(17)

where k stands for the Boltzmann’s constant and $\nu$ for the frequency of the (moving) ion vibrating in the direction of the potential barrier. By accepting the approximation [23] that $\nu ={\nu }_{TO}(k\to 0)$, where ${\nu }_{TO}$ is the frequency at very long wavelengths of the transverse optical mode, and using in Equation (17) the values ${\nu }_{TO}=0.92\times {10}^{13} {\mathrm{s}}^{-1}$ [24] and ${\tau }_0^{-1}=1.9\times {10}^{14} {\mathrm{s}}^{-1}$, we find $s_{bv}^m=2.3 \mathrm{k}$.

Since both parameters $h_{bv}^m=0.66$ eV and $s_{bv}^m=2.3$ k are now known, the $g_{bv}^m$ value can be obtained by means of $g_{bv}^m=h_{bv}^m-\mathrm{T}{s}_{bv}^m$ at any temperature. Hence, at T = T${}_{max}$ = 211 K, we find $g_{bv}^m$ = 0.62 eV. The aforementioned value $s_{bv}^m=2.3$ k is comparable to the one ($s^m=2.3$ k) of the free cation vacancy motion identified in the measurements in Ref. [25], which also resulted in $h^m$ = 0.75 eV leading to $g^m$ = 0.69 eV at T = 300 K for the cation vacancy motion.

As far as the Schottky defect formation process is concerned, Lallemand [25] found $h^f$ = 2.37 eV and $s^f$ = 8.9 k, which result in $g^f$ = 2.14 eV at T = 300 K. The experimental results for the defect volumes are $v^f$ = (12.6–13.2) cm${}^3$/mole and $v^m$ = (3.2–4.0) cm${}^3$/mole for the cation vacancy migration process [25], and hence the corresponding $\frac{v^i}{g^i}$ values are

$$\frac{v^f}{g^f}=(5.9-6.2)\times {10}^{-2} {\mathrm{GPa}}^{-1}$$

(18)

and

$$\frac{v^m}{g^m}=(4.6-5.8)\times {10}^{-2} {\mathrm{GPa}}^{-1}$$

(19)

for these two defect processes, respectively.

Assuming that the defect migration volume $v_{bv}^m$ for the bound cation vacancy motion is approximately equal [12] to that of the free cation vacancy motion, we find

$$\frac{v_{bv}^m}{g_{bv}^m}=(5.2-6.5)\times {10}^{-2} {\mathrm{GPa}}^{-1}$$

(20)

According to Equation (15), the ratio $\frac{v^i}{g^i}$ should be equal to $\frac{1}{B}\left({\left.\frac{dB}{dP}\right|}_T-1\right)$ irrespective of the defect processes. By considering the values B = 62.3 GPa and ${\left.\frac{dB}{dP}\right|}_T$ = 5.01 reported in the third column of Table 1 of Dewaele [20], as mentioned, we find

$$\frac{1}{B}\left({\left.\frac{dB}{dP}\right|}_T-1\right)=6.4\times {10}^{-2} {\mathrm{GPa}}^{-1}$$

(21)

By the same token, the values B = 66.2 GPa and ${\left.\frac{dB}{dP}\right|}_T$ = 4.6 obtained by Dong et al. [19] lead to

$$\frac{1}{B}\left({\left.\frac{dB}{dP}\right|}_T-1\right)=5.4\times {10}^{-2} {\mathrm{GPa}}^{-1}$$

(22)

In order to visualize the agreement between the experimental $\frac{v^i}{g^i}$ values and the bulk quantity $\frac{1}{B}\left({\left.\frac{dB}{dP}\right|}_T-1\right)$ predicted by the cB$\mathrm{\Omega }$ model, the $v^i$ values are plotted in Figure 2 versus the $g^i$ values. In the same figure, we also plot the straight lines corresponding to the predicted values $\frac{1}{B}\left({\left.\frac{dB}{dP}\right|}_T-1\right)$ of the ratio $\frac{v^i}{g^i}$ from Equation (15) deduced either from the above-mentioned measurements by Dewaele [20] (blue straight line) or by Dong et al. [19] (magenta straight line).

Figure 2 shows that, for all three defect processes, the experimental results of $\frac{v^i}{g^i}$ agree with those deduced from the cB$\mathrm{\Omega }$ model. A more detailed inspection of Figure 2 uncovers that the experimental results of $\frac{v^i}{g^i}$ lie closer to the blue straight line, which remarkably corresponds to the latest equation of state that appeared in Ref. [20].

In short, the experimental results of various defect processes in LiF confirm the cB$\mathrm{\Omega }$ model upon using the most recent equation of state [20].

4. Applications of the cB$\mathrm{\Omega }$ Model to Nuclear Fuels and Materials for Microelectronic Devices

The strength of the cB$\mathrm{\Omega }$ model to describe the defect processes of nuclear and electronic materials has been recently highlighted considering, for example, nuclear fuel [26,27,28] and materials for microelectronics devices [29,30,31]. Regarding nuclear fuel, there is considerable work on UO${}_2$; however, the advent of proliferation resistant fuel such as ThO${}_2$ and novel fuels (mixed oxide) requires a systematic understanding of their material properties in a wide range of conditions (pressure and temperature). Cooper et al. [27] investigated oxygen self-diffusion in UO${}_2$ using molecular dynamics in conjunction with the cB$\mathrm{\Omega }$ model. This theoretical approach allowed the study to span over a considerable range of temperatures (300–1900 K) and pressures (0–10 GPa), which would be difficult and expensive to consider experimentally. Importantly, it was highlighted that there is a considerable reduction in oxygen self-diffusion (increase in the activation energy of oxygen self-diffusion) as a function of the increasing pressure.

The above study by Chroneos and coworkers in nuclear fuels was further elaborated in Ref. [32], which revealed the very interesting result that thermal expansivity and compressibility of ${\upsilon }^{act}$ in the oxygen self-diffusion in UO${}_2$ exceed the corresponding bulk quantities almost by one order of magnitude, see Figure 3 and Figure 4 for the ratios ${\beta }^{act}/\beta$ and ${\kappa }^{act}/\kappa$ versus temperature for several pressures. More precisely, these results have been obtained as follows:

Cooper et al. [27] derived both the expansivity and B using molecular dynamics (MD) calculations and found the following expression for B(T,P):

$$B(T,P)=a_0+b_{0}T+c_0{T}^2+d_{0}P+e_0{P}^2+f_{0}PT$$

(23)

where $a_0$ = 218.0 GPa is the bulk modulus at T = 0, P = 0, $b_0 = -4.33\times {10}^{-2}$ GPa K${}^{-1}$, $c_0 = -1.846\times {10}^{-6}$ GPa K${}^{-2}$, $d_0$ = 5.864 is the pressure derivative of B, i.e., $\frac{dB}{dP}$, at T = 0, P = 0, $e_0 = -1.387\times {10}^{-1}$ GPa${}^{-1}$ and $f_0=1.301\times {10}^{-3}$ K${}^{-1}$ (see Table 1 of Cooper et al. [27]). By means of Equation (23), the use of Equation (10) leads to the ratio ${\kappa }^{act}/\kappa$ at various temperatures for various pressures depicted in Figure 3. This figure reveals that ${\kappa }^{act}/\kappa$ exceeds unity, as expected from Equation (10) because in its last term we have [12]: $\frac{d^{2}B}{d{P}^2}<0$ in the numerator and $\frac{dB}{dP}{|}_T>1$ in the denominator. We find that ${\kappa }^{act}$ exceeds $\kappa$ by almost an order of magnitude or so, which is important because it is usually assumed [33,34,35] that the activation volume changes with pressure similarly with the crystal volume, i.e., ${\kappa }^{act}=\kappa$ [36]. Furthermore, we also calculate ${\beta }^{act}/\beta$ at various temperatures and pressures by applying Equation (11), and the results are depicted in Figure 4, where $\beta$ has been calculated from the relation:

$$\mathrm{\Omega }={\mathrm{\Omega }}_0+a_{1}P+a_2{P}^2+{\beta }_{1}T+{\beta }_2{T}^2+{\gamma }_{0}PT$$

(24)

The latter relation is similar to that of Equation (23) and results from a polynomial fit to the data presented in Figure 1 of Ref. [27] yielding ${\mathrm{\Omega }}_0=13.513$ Å${}^3$, $a_1=-6.387\times {10}^{-2}$ Å${}^3$ GPa${}^{-1}$, $a_2=1.267\times {10}^{-3}$ Å${}^3$ GPa${}^{-2}$, ${\beta }_1=3.881\times {10}^{-4}$ Å${}^3$ K${}^{-1}$, ${\beta }_2=4.6\times {10}^{-8}$ Å${}^3$ K${}^{-2}$, and ${\gamma }_0=-1.662\times {10}^{-5}$ Å${}^3$ GPa${}^{-1}$ K${}^{-1}$; it gives B values compatible within ±5% to the ones obtained when using Equation (23) for $P<6$ GPa and $T\le 1360$ K. Figure 4 shows that ${\beta }^{act}/\beta$ exceeds unity (by almost one order of magnitude) as expected from Equation (11) if we consider that the numerator of the last term is positive since the value of $\frac{dB}{dP}{|}_T$ increases upon increasing the temperature due to anharmonic effects [12]. Concerning Figure 3 and Figure 4, the following comments are now in order: These two figures have been plotted in the temperature 300–1500 K, which somewhat exceeds (i.e., by around 140 degrees) the higher temperature limit $T\le 1360$ K for the elastic and expansivity data mentioned above when using Equations (23) and (24). This does not change the importance of the results obtained here, i.e., ${\kappa }^{act}/\kappa$ as well as ${\beta }^f/\beta$ exceed unity by almost one order of magnitude in view of the following: The curvature of the Arrhenius plots (with one operation mechanism) is quantified by the quantity: $\frac{d{h}^{act}}{dT}{\mid }_P=\frac{T{\upsilon }^{act}{\beta }^2}{\kappa }(2\frac{{\beta }^{act}}{\beta }-\frac{{\kappa }^{act}}{\kappa })+\mathrm{a}\phantom{\mathrm{x}}\mathrm{small}\phantom{\mathrm{x}}\mathrm{term},$ see p. 388 of Ref. [12]. We find that the quantity between parentheses changes only slightly from the value 5.3 to 5.2 for the aforementioned two cases $T=1360$ K and 1500 K, respectively, which means that the extent of the Arrhenius curvature for oxygen self-diffusion plots in UO${}_2$ studied by Cooper et al. [27] does not seriously change upon investigating the above two temperature ranges.

Considering materials for electronics, the aggressive scaling for decades has led to characteristic dimensions of a few nanometers. This in turn requires the understanding of intrinsic defect properties and dopant-defect interaction in both established materials (silicon) [30,31] and emerging materials (germanium, silicon germanium, and III-V alloys) [29,31]. For example, Kube et al. [37] determined a curved Arrhenius plot on self-diffusion in Si and proposed that it is associated with two diffusion mechanisms (vacancies and self-interstitials). On the contrary, on the basis of the cB$\mathrm{\Omega }$ model, Saltas et al. [30] proposed that the experimental results can also be comprehended by a single diffusion mechanism if the nonlinear anharmonic behaviour of the isothermal bulk modulus is considered. Interestingly, the analysis on the ground of the cB$\mathrm{\Omega }$ model reveals that the activation enthalpy of self-diffusion varies nonlinearly from 3.6 eV to 4.9 eV in good agreement with the experimental data by Kube et al. [37]. Finally, there are numerous studies on germanium, silicon germanium, and III-V alloys and, in these systems, the cB$\mathrm{\Omega }$ model can provide complementary information to the limited experimental data (refer to Ref. [38] and references therein).

5. Application of the cB$\mathrm{\Omega }$ Model to High T${}_{\mathit{c}}$ Superconductors

In high T${}_c$ superconductors, like $YB{a}_{2}C{u}_3{O}_{7-\delta }$, high critical current densities can be achieved by the presence of high-density defects which will provide suitable pinning centers for the magnetic flux lines (e.g., see Ref. [39] and references therein). In such a flux-line pinning, point defects play an important role in $YB{a}_{2}C{u}_3{O}_{7-\delta }$ in which the superconductivity coherence length is of order of tens of Å. Doping $YB{a}_{2}C{u}_3{O}_{7-\delta }$ with alkaline earth elements improve bulk and grain-boundary transport [40] and in addition was found by Su et al. [41] to lead to a formation energy for a Schottky defect, which interestingly is in accordance with the expectations of the cB$\mathrm{\Omega }$ model (e.g., see [8,10,26,38,42]). This model, however, does not seem to describe the case when replacing Y in the doping of $YB{a}_{2}C{u}_3{O}_{7-\delta }$ by rare earth elements (which is strikingly reminiscent, of an earlier finding [43] that in alkali halides doped with alkaline earth divalent cations, the $h^m$ value for the reorientation process of the electric dipoles formed by a divalent cation vacancy increases upon increasing the ionic radius of the dopants but not for rare earth divalent dopant). In a work which just appeared [44], the influence of high hydrostatic pressure up to 11 kbar on the conductivity in the basal ab-plane of medium-doped with praseodymium ($x\approx 0.23$) single-crystal $Y_{1-x}P{r}_{x}B{a}_{2}C{u}_3{O}_{7-\delta }$ samples was investigated and the following interesting fact emerged: in contrast to the pure $YB{a}_{2}C{u}_3{O}_{7-\delta }$ samples with the optimal oxygen content, the application of high pressure leads to the formation of phase separation in the basal plane of $Y_{0.77}P{r}_{0.23}B{a}_{2}C{u}_3{O}_{7-\delta }$ single crystals. Possible mechanisms of the effect of Pr doping and high pressure on the two-step resistive transition to the superconducting state were discussed, and it was determined that, in the normal state, the conductivity is metallic limited by phonon scattering (Bloch-Grüneisen regime) and defects. It has been concluded that doping with Pr leads to the appearance of an anisotropic distribution of defects and flat macroscopic superconducting regions with different T${}_c$ are formed. A concentration increase of Pr leads to a decrease in the localization length of charge carriers, which indicates a decrease in the spatial propagation of the wave function of charge carriers. In the normal state, the main contribution to the resistivity of the $Y_{0.77}P{r}_{0.23}B{a}_{2}C{u}_3{O}_{7-d}$ single crystal is made by intraband scattering ($\propto T^5$) and residual resistivity. Concerning the hydrostatic pressure, it leads to anisotropy decrease of linear compressibility along with a decrease in phonon resistivity and residual resistivity, while Debye temperature and coherence length remain constant.

6. Applications of Usefulness in Geophysics and Seismology

We start with the following three cases: self-diffusion in diamond, self-diffusion in iron and the thermodynamic estimation of the compressibility of ferropericlase under high pressure.

In earth sciences, diamond diffusion properties are important because natural diamonds and their mineral inclusions provide information about the geochemical character and geotherm of the ancient continental lithosphere e.g., see Ref. [45] and references therein. Specifically, the spatial distribution of carbon and nitrogen isotopes in diamond provides information on mantle residence time [46]. Diamonds with long residence at high temperature will gradually lose their initial zoning patterns due to diffusion; thus, the diffusion data can constrain the maximum possible age of diamonds. Rare diamonds that originate from the mantle transition zone should have developed [46] a length scale of ≈1 mm of isotopic zoning over the age of the Earth $4.5\times {10}^9$ years. A calculation of point defect parameters in diamond based on the cB$\mathrm{\Omega }$ model has been initially made by Varotsos [39] in 2007, which was followed in 2012 by a calculation of the self-diffusion coefficients by Zhang and Wu [47]. These authors, in light of more recent pressure–volume–temperature relationships, showed that the self-diffusion coefficient of diamond as a function of temperature and pressure can be reproduced by means of the cB$\mathrm{\Omega }$ model using bulk elastic and expansivity data. In particular, Zhang and Wu stated that they found an excellent agreement between the calculated self-diffusion coefficients ($D_{calc}$) over a range of temperatures (1000–2800 K) and pressures (0–80 GPa), see their Figure 1. Specifically, they found from the slope and intercept ($D_o$) of $ln{D}_{calc}$ vs. $1/T$ the values $h^{act}=6.71\pm 0.20$ eV and $D_o=3.{30}_{-0.21}^{+0.19}\times {10}^{-5} {\mathrm{m}}^2/\mathrm{s}$ which are comparable to the experimental derived values $6.8\pm 1.6$ eV and $4.1_{-3.5}^{+26.3}\times {10}^{-5} {\mathrm{m}}^2/\mathrm{s}$, respectively.

Concerning iron, it is considered as the major component of the Earth’s core [48] and since it is known that the diffusion is related strongly to the creep, grain boundary motion, and solid phase transition in a metal and alloy, a knowledge of the transport properties of iron and its alloys is important for understanding the evolution and current dynamics of the cores of many of the planets, including Earth. In 2014, Zhang [49] used the available $P-V-T$ equation of state of iron and successfully reproduced the temperature and pressure dependence of self-diffusion coefficients in iron polymorphs ($\alpha ,\delta ,\gamma$ and $ϵ$ phases) in terms of the bulk elastic and expansivity data by means of the cB$\mathrm{\Omega }$ model. It was found that the calculated diffusion parameters, such as self-diffusion coefficient, activation energy and activation volume, over a broad temperature range (500–2500 K) and pressure range (0–100 GPa), compare favorably well with the experimental ones when the uncertainties are considered. Interestingly, this approach has been satisfactorily applied to the case of $ϵ$-Fe under the corresponding $P-T$ (up to 100 GPa and 2500 K) conditions of the lower mantle; such data should have important implications for the physical and chemical processes of deep planetary interiors [50,51].

We now turn to ferropericlase (Mg,Fe)O, which is widely believed to be a major phase in the lower mantle of the Earth, coexisting with Mg-rich silicate perovskite, (Mg,Fe,Al)SiO${}_3$. Since the elastic properties of (Mg${}_{1-x}$Fe${}_x$)O ferropericlase are essential to analyze seismic data and to constrain its chemical composition in the lower mantle, Zhang and Xu [52] proceeded to the estimation of the bulk moduli of (Mg${}_{1-x}$Fe${}_x$)O ferropericlase as a function of component x in terms of the elastic data of the end members. Along these lines, they followed the method proposed by Varotsos and Alexopoulos [12] to explain how the compressibility $\kappa (=1/B)$ of a solid solution A${}_{1-x}$B${}_x$ is interrelated with the compressibilities of two end members A and B upon using the cB$\mathrm{\Omega }$ model. Zhang and Xu [52] stated that they found an excellent consistency between the calculated B values and the experimental ones in ferropericlase. Hence, they concluded that this method based on the cB$\mathrm{\Omega }$ model exhibits credibility of describing compresibility in terms of the bulk moduli of the two end members, i.e., MgO and FeO; thus, it constitutes an interesting method to determine the compressibility of Earth materials with complicated structure using the measured compressibilities of their end members.

Finally, we point out a mechanism that has been proposed in Ref. [53], which achieves in bridging the water solubility and ion diffusivity in the Earth’s mantle silicates by using the cB$\mathrm{\Omega }$ model. This is important because, since numerous studies have demonstrated that even tens of ppm of water by weight also exerts a considerable effect on the transport properties of mantle minerals such as diffusion [54,55,56,57], rheology (Mei and Kohlstedt [58]; Xu et al. [59]), and electrical conductivity [60,61,62], an understanding of water dissolution in nominally anhydrous minerals at the atomic level is fundamental for the development of dynamic models of the Earth’s mantle. This is crucial because several theoretical and experimental studies have revealed that water can be incorporated into silicate minerals by a variety of potential defect mechanisms, such as metal site vacancy, Si site vacancy, and O ions [63,64,65,66,67,68], which makes the role of the cB$\mathrm{\Omega }$ model important in this case.

Beyond the aforementioned three examples, several other applications of the cB$\mathrm{\Omega }$ model to Geophysics have been published, e.g., calculation of diffusion parameters of He in olivine [69], calculation of oxygen self-diffusion coefficients in minerals [70], calculation of diffusion parameters of Si in silicates [71] or of Si self-diffusion in stishovite [72]; for more such applications, see Ref. [38].

7. Application of the cB$\mathrm{\Omega }$ Model to a Problem of Major Technological Interest

Si-based anodes can drastically improve the energy density of current in Li-ion batteries, but their actual use is still limited by a huge volume change (up to 300%) during charging and discharging [73]. This is due to complex structural transformations that Si experiences during Li insertion and removal [74]. It is therefore of paramount importance to investigate the kinetics of Li in both bulk Si and various lithiated phases LiSi in order to understand the structural changes.

An approach to obtain a non-classical description of diffusion based on the transition state theory (TST) has been proposed [75,76], which is often used with the harmonic approximation for computing the phonon frequencies; thus, it is important to estimate the effect of anharmonicity in the system that might become significant at a high temperature. Such a methodology based on TST was recently followed by Kumar et al. [77] to study Li interstitial diffusion in bulk Si, as well as Li diffusion in the crystalline lithiated phase LiSi. They also took into account nuclear quantum approximation and anharmonic effects. Concerning the latter, Kumar et al. [77] pointed out that, in the literature on anharmonicity in diffusion, mainly two types of effects are discussed: (1) those due to the change in the lattice parameter owing to thermal expansion at a high temperature [78,79] and (2) those due to the phonon–phonon interactions [80,81]. The first effect, i.e., the volume change, is often studied through the quasiharmonic approximations that provide a first order anharmonic correction [82]. However, it can underestimate the phonon frequency change when the phonon–phonon interactions contribute significantly in phonon softening at high temperatures [80]. In order to address both effects, Kumar et al. [77] considered the cB$\mathrm{\Omega }$ model. They finally concluded that, as compared to the classical approach, the addition of quantum effects for Li diffusion in Si lowers the value of the room temperature diffusion coefficient by 33% because of a higher effective energy barrier. The cB$\mathrm{\Omega }$ model and the quasi-harmonic approximation suggest an increase in the diffusion coefficient by 60% in the case of bulk Si by including anharmonic effects. In the case of LiSi, Kumar et al.’s calculations are in better agreement with the experimental data than previous ab initio molecular dynamics studies conducted at higher temperatures. They also obtained an effective diffusion barrier of 0.27 eV ± 0.01 eV with a pre-factor of 2.8 $\pm 1.1\times {10}^{-5} {\mathrm{cm}}^2/\mathrm{s}$ for Li diffusion in LiSi.

8. Search for Composition of Better Target Properties in Cu-Co-Si Alloys via Machine Learning

Zhao et al., in a very recent work [83], limiting the discovery of alloys with high hardness and good electrical conductivity (EC), propose the combination of a Gaussian process regression-based machine learning (ML) surrogate model and the phase diagram to search for compositions with better target properties in Cu-Co-Si alloys. Since Ag is more expensive, Cu finds wide use in household and industrial products in view of its excellent conductivity. Other alloying elements are introduced to pure Cu to improve its strength [84,85]. The introduction of alloying elements brings strengthening of Cu alloys and in particular Si is one of the common alloying elements in Cu alloys because it produces a strong solution strengthening effect and decreases the surface tension as well as improves the fluidity of Cu melt [86]. However, the equilibrium solubility of Si in Cu lattice is as high as 4.6 wt% at an ageing temperature of 450 ${}^{\circ }$C [87], which contributes to strong electron scattering and a decrease in EC. Co can react with Si to form CoSi and Co${}_2$Si compounds. In addition, the equilibrium solubility of Co in the Cu matrix at 450 ${}^{\circ }$C is close to 0 wt% [88]. Therefore, the addition of Co has the potential to improve the strength and conductivity of the Cu-Si alloy simultaneously.

As for ageing strengthened Cu alloys, the precipitation phase formed by alloying elements affects the conductivity and hardness. To confirm the influence of the compositions and microstructures on properties in Cu-Co-Si alloys, their ML model automatically designed 12,000 alloys and predicted their properties. Combined with the phase diagram and microstructure characterization, they concluded that the EC and hardness of Cu-Co-Si alloys are related to the moderate Co${}_2$Si and CoSi precipitation.

For metals with vacancy diffusion as the main diffusion mechanism, atomic diffusion requires diffusion activation Gibbs energy to provide power. According to the cB$\mathrm{\Omega }$ model [6,11], the defect Gibbs energy is proportional to the isothermal bulk modulus, thus the diffusion coefficient D was obtained from Equation (13) (see also Section 2).

The formation of precipitation driven by atom diffusion is the core reason for the high strength and high conductivity of ageing strengthened copper alloys. Therefore, Zhao et al. [83] inferred that the isothermal bulk modulus and compressibility are key factors to predict hardness and electrical conductivity. This is why the cB$\mathrm{\Omega }$ model plays an important role in the aforementioned research for composition of better target properties in Cu-Co-Si alloys.

9. Recent Studies on the Foundation of the Physical Origin of the cB$\mathrm{\Omega }$ Model

During the last two years, two studies appeared which shed more light on the foundation of the physical origin of the cB$\mathrm{\Omega }$ model. They are summarised below:

9.1. Foundation of cB$\mathrm{\Omega }$ Model Studying the Interfacial Mobility Gradient in Thin Polymer Films Using Molecular Dynamics Simulations on a Layer-by-Layer Basis

In Ref. [89], Zhang et al. examined the mobility gradient in the interfacial region of substrate-supported polymer films using molecular dynamics simulations and interpreted these gradients within the string model of glass-formation. After analyzing the mobility gradient on a layer-by-layer basis, they found that no large gradients in the extent of collective motion exist in these simulated films. In other words, the novelty of this work by Zhang et al. [89] is that they applied the string model on a layer-by-layer basis in the supported polymer film, which allowed the quantification of the activation parameters $h^{act}(z)$ and $s^{act}(z)$ on a layer-by-layer basis. They concluded that the interfacial mobility gradient derives mainly from a gradient in the high-temperature activation enthalpy $h^{act}$ and entropy $s^{act}$ as a function of depth z. In order to gain more insight into the interfacial mobility gradient, they also examined various material properties suggested previously to influence $h^{act}$ in condensed materials, including density, potential, cohesive energy density, and a local measure of stiffness or ${\langle u^2(z)\rangle }^{-3/2}$, where $\langle u^2(z)\rangle$ is the average mean squared particle displacement at a caging time on the order of a ps. They also finally found that changes in local stiffness best correlate with changes in $h^{act}(z)$ and that $s^{act}(z)$ contributes significantly to the interfacial mobility gradient, so it must not be neglected. This constitutes a direct confirmation of the cB$\mathrm{\Omega }$ model because it clearly demonstrates that $h^{act}$ and $s^{act}$ best correlate with changes of the bulk modulus B as follows (note that, for example, Equation (6) can be alternatively written as $h^{act}=c^i\mathrm{\Omega }{B}_0^{SL}$, where $B_0^{SL}$ is the intercept of the linear part of the graph $B=f(T)$ with the vertical axis since it is found experimentally that B decreases almost linearly upon increasing temperature, e.g., see Equations (5)–(12) in Ref. [4]): In Ref. [89], it was clearly stated that earlier studies [90,91] have shown a quantitative relation between the bulk modulus and Debye–Waller (DW) factor $\langle u^2(z)\rangle$, suggesting that $\langle u^2(z)\rangle$ at a local scale might provide a local measure of material stiffness. The DW factor $\langle u^2(z)\rangle$ is defined by the mean square particle displacement at a caging time on the order of a ps. Physically, one can expect that the amplitude of molecular vibrations is affected by the local stiffness of the material, providing some intuition for the possible relationship of $\langle u^2(z)\rangle$ and the bulk modulus. At a more basic level, $\langle u^2(z)\rangle$ as a measure of molecular rigidity derives from the harmonic oscillator model in which the harmonic oscillator force constant $k_h$ equals $kT/\langle u^2(z)\rangle$ [92]. In short, Ref. [89] of all the quantities considered, ${\langle u_A^2(z)\rangle }^{-3/2}$ (which they argued relates to bulk modulus) best tracks the variation of $h^{act}$ with film depth.

9.2. Analytical Method for the Foundation of the cB$\mathrm{\Omega }$ Model Using Mie–Lennard–Jones Interatomic Potential

In Ref. [93], Magomedov proposed an analytical method for calculating the parameters of vacancy formation and self-diffusion of atoms in a single-component crystal. The method is based on the four-parameter pairwise Mie–Lennard–Jones interatomic interaction potential. The method allows calculation of all activation process parameters, that is, Gibbs energy, enthalpy, entropy, and volume, for both the vacancy formation process and the self-diffusion process. As an example, Magomedov presented the case for gold calculated from T = 10–1330 K along two isobars, P = 0 and 24 GPa. In the high-temperature region, the vacancy formation and the self-diffusion coefficient pass into classical Arrhenius dependencies, with a weakly temperature-dependent enthalpy and a positive value for the activation process entropy, which is a behavior strikingly similar with that exhibited by the cB$\mathrm{\Omega }$ model, as it was published since 1977 [2]. In short, Magomedov’s analytical method leads to results that are fully compatible with the cB$\mathrm{\Omega }$ model.

10. Summary and Conclusions

In the 1980s, a strict foundation of thermodynamics of point defects in solids appeared and revealed the existence of two families of defect parameters pointing to the key role of isothermal bulk modulus. Thus, the cB$\mathrm{\Omega }$ model naturally emerged as a tentative interrelation between point defect parameters and bulk properties. Here, the following results have been obtained:

(i) During the last decade, in a multitude and variety of fields, several successful applications of the cB$\mathrm{\Omega }$ model appeared including high pressure diamond anvil measurements, high T${}_c$ superconductors, nuclear fuels, and materials for micro-electronics devices, applications of usefulness in Geophysics and Seismology, a problem of major technological interest and search for compositions of better target properties in Cu-Co-Si alloys via machine learning.

(ii) During the last few years, the foundation of the physical origin of the cB$\mathrm{\Omega }$ model has been further strengthened by means of two independent studies.

Hence, in view of the above recent results, we conclude that such an interconnection between the point defect parameters and bulk properties in solids does exist.