Phonon localization in binary alloys with diagonal and off-diagonal disorder: A cluster Green's function approach

TL;DR

This paper introduces a new Green's function-based method to study phonon localization in disordered alloys, revealing how mass and force-constant disorder influence phonon spectra and localization, with implications for thermoelectric materials.

Contribution

The paper develops and validates a novel TMDCA-based computational approach for analyzing phonon spectra in disordered alloys, incorporating both mass and force-constant disorder effects.

Findings

Boson peak emerges at low soft particle concentrations.

States at the boson peak crossover from localized to extended with increasing soft particle concentration.

Vacancies strongly localize phonons even at low concentrations.

Abstract

We report the development and application of a new method for carrying out computational investigations of the effects of mass and force-constant (FC) disorder on phonon spectra. The method is based on the recently developed typical medium dynamical cluster approach (TMDCA), which is a Green's function approach. Excellent quantitative agreement with previous exact diagonalization results establishes the veracity of the method. Application of the method to a model system of binary mass and FC-disordered system leads to several findings. A narrow resonance, significantly below the van Hove singularity, that has been termed as the boson peak, is seen to emerge for low soft particle concentrations. We show, using the typical phonon spectrum, that the states constituting the boson peak cross over from being completely localized to being extended as a function of increasing soft particle concentration. In general, an interplay of mass and FC disorder is found to be cooperative in nature, enhancing phonon localization over all frequencies. However, for certain range of frequencies, and depending on material parameters, FC disorder can delocalize the states that were localized by mass disorder, and vice-versa. Modeling vacancies as weakly bonded sites with vanishing mass, we find that vacancies, even at very low concentrations, are extremely effective in localizing phonons. Thus, inducing vacancies is proposed as a promising route for efficient thermoelectrics. Finally, we use model parameters corresponding to the alloy system, Ni1-xPtx, and show that mass disorder alone is insufficient to explain the pseudogap in the phonon spectrum; the concomitant presence of FC disorder is necessary.

Equations

$$H=\sum _{l\alpha }\frac{p_{\alpha }^2(l)}{2M(l)}+\frac{1}{2}\sum _{l{l}^{\prime },\alpha \beta }{\Phi }^{\alpha \beta }(l,l^{\prime })u_{\alpha }(l)u_{\beta }(l^{\prime }),$$

$$\begin{array}{rl}{x}_l=1,y_l=0,& \text{if l∈A}\\ x_l=0,y_l=1,& \text{if l∈B}\end{array}$$

$$\begin{array}{r}{x}_l{y}_l=0,x_l^2=x_l\\ ⟨x_l⟩=c_{\mathrm{A}},⟨y_l⟩=c_{\mathrm{B}}\end{array}$$

$$M(l)=\{\begin{array}{l}{x}_{l}M(l)x_l=M_{\mathrm{A}}\\ y_{l}M(l)y_l=M_{\mathrm{B}}\\ x_{l}M(l)y_l=M_{\mathrm{AB}}=0\\ y_{l}M(l)x_l=M_{\mathrm{BA}}=0\end{array}$$

$${\left(\hat{V}\right)}_{l{l}^{\prime }}=\left(1-M(l)/M_0\right){\delta }_{l,l^{\prime }},$$

$$P(V_l)=c_{\mathrm{A}}\delta (V_l-V_{\mathrm{A}})+c_{\mathrm{B}}\delta (V_l-V_{\mathrm{B}}),$$

$$\Phi (l,l^{\prime })=\{\begin{array}{l}{x}_l{\Phi }^{\alpha \beta }(l,l^{\prime })x_{l^{\prime }}={\Phi }^{\alpha \beta ,\mathrm{AA}}(l,l^{\prime })\\ y_l{\Phi }^{\alpha \beta }(l,l^{\prime })y_{l^{\prime }}={\Phi }^{\alpha \beta ,\mathrm{BB}}(l,l^{\prime })\\ x_l{\Phi }^{\alpha \beta ,}(l,l^{\prime })y_{l^{\prime }}={\Phi }^{\alpha \beta ,\mathrm{AB}}(l,l^{\prime })\\ y_l{\Phi }^{\alpha \beta }(l,l^{\prime })x_{l^{\prime }}={\Phi }^{\alpha \beta ,\mathrm{BA}}(l,l^{\prime })\end{array}$$

$${\Phi }^{\alpha \beta }(l,l^{\prime })={\delta }_{\alpha \beta }({\Phi }_D{\delta }_{l,l^{\prime }}+{\Phi }_{nn}{\delta }_{{\mathbf{R}}_{l^{\prime }},{\mathbf{R}}_l+\vec{\delta }}),$$

$$M(l){\omega }^{2}D(l,l^{\prime },\omega )={\delta }_{l{l}^{\prime }}+\sum _{l^{\prime \prime }{l}^{\prime }}\Phi (l,l^{\prime \prime })D(l^{\prime \prime },l^{\prime },\omega )$$

$$D_{\mathrm{AA}}(l,l^{\prime })$$

$$+d_{\mathrm{A}}(l)\sum _{l^{\prime \prime }\ne l}{\Phi }^{\mathrm{AB}}(l,l^{\prime \prime })D_{\mathrm{BA}}(l^{\prime \prime },l^{\prime })$$

$$D_{\mathrm{AB}}(l,l^{\prime })$$

$$+d_{\mathrm{A}}(l)\sum _{l^{\prime \prime }\ne l}{\Phi }^{\mathrm{AB}}(l,l^{\prime \prime })D_{\mathrm{BB}}(l^{\prime \prime },l^{\prime })$$

$$D_{\mathrm{BA}}(l,l^{\prime })$$

$$+d_B(l)\sum _{l^{\prime \prime }\ne l}{\Phi }^{\mathrm{BB}}(l,l^{\prime \prime })D_{\mathrm{BA}}(l^{\prime \prime },l^{\prime })$$

$$D_{\mathrm{BB}}(l,l^{\prime })$$

$$+d_{\mathrm{B}}(l)\sum _{l^{\prime \prime }\ne l}{\Phi }^{\mathrm{BB}}(l,l^{\prime \prime })D_{\mathrm{BB}}(l^{\prime \prime },l^{\prime })$$

$$x_{l}D(l,l^{\prime })x_{l^{\prime }}=D_{\mathrm{AA}}(l,l^{\prime })$$

$$x_{l}D(l,l^{\prime })y_{l^{\prime }}=D_{\mathrm{AB}}(l,l^{\prime })$$

$$y_{l}D(l,l^{\prime })x_{l^{\prime }}=D_{\mathrm{BA}}(l,l^{\prime })$$

$$\displaystyle y_{l}D(l,l^{\prime})y_{l^{\prime}}=D_{\rm BB}(l,l^{\prime})\$$

$$x_{l}d(l)x_l$$

$$y_{l}d(l)y_l$$

$$\begin{gathered} {\left(\begin{array}{cc}{D}_{\mathrm{AA}}& D_{\mathrm{AB}}\\ D_{\mathrm{BA}}& D_{\mathrm{BB}}\end{array}\right)}_{l{l}^{\prime }}={\left(\begin{array}{cc}{d}_{\mathrm{A}}& 0\\ 0& d_{\mathrm{B}}\end{array}\right)}_l{\delta }_{l{l}^{\prime }} \\ + \end{gathered}$$

$$\begin{gathered} {\left(\begin{array}{cc}{d}_{\mathrm{A}}& 0\\ 0& d_{\mathrm{B}}\end{array}\right)}_l \\ \sum _{l^{\prime \prime }\ne l}{\left(\begin{array}{cc}{\Phi }^{\mathrm{AA}}& {\Phi }^{\mathrm{AB}}\\ {\Phi }^{\mathrm{BA}}& {\Phi }^{\mathrm{BB}}\end{array}\right)}_{l{l}^{\prime \prime }} \\ {\left(\begin{array}{cc}{D}_{\mathrm{AA}}& D_{\mathrm{AB}}\\ D_{\mathrm{BA}}& D_{\mathrm{BB}}\end{array}\right)}_{l^{\prime \prime }{l}^{\prime }} \end{gathered}$$

$$\underline{D}=\underline{d}+\underline{d}\times \underline{\Phi }\times \underline{D}$$

$$\underline{\Delta (\mathbf{K},\omega )}=\left(\begin{array}{cc}{\Delta }^{\mathrm{AA}}(\mathbf{K},\omega )& {\Delta }^{\mathrm{AB}}(\mathbf{K},\omega )\\ {\Delta }^{\mathrm{BA}}(\mathbf{K},\omega )& {\Delta }^{\mathrm{BB}}(\mathbf{K},\omega )\end{array}\right)$$

$$\Phi^{\prime}(l,l^{\prime})=\begin{cases}\sum_{\mathbf{K}}{\big{(}\underline{\omega^{2}_{{\mathbf{K}}}}\big{)}^{\rm AA}}e^{i\mathbf{K}\cdot(\vec{R}_{l}-\vec{R}_{l^{\prime}})},&\text{if $l\in\rm A$, $l^{\prime}\in\rm A$ }\\ \sum_{\mathbf{K}}{\big{(}\underline{\omega^{2}_{{\mathbf{K}}}}\big{)}^{\rm AB}}e^{i\mathbf{K}\cdot(\vec{R}_{l}-\vec{R}_{l^{\prime}})},&\text{if $l\in\rm A$, $l^{\prime}\in\rm B$ }.\\ \sum_{\mathbf{K}}{\big{(}\underline{\omega^{2}_{{\mathbf{K}}}}\big{)}^{\rm BA}}e^{i\mathbf{K}\cdot(\vec{R}_{l}-\vec{R}_{l^{\prime}})},&\text{if $l\in\rm B$, $l^{\prime}\in\rm A$ }.\\ \sum_{\mathbf{K}}{\big{(}\underline{\omega^{2}_{{\mathbf{K}}}}\big{)}^{\rm BB}}e^{i\mathbf{K}\cdot(\vec{R}_{l}-\vec{R}_{l^{\prime}})},&\text{if $l\in\rm B$, $l^{\prime}\in\rm B$ }.\end{cases}$$

$$\underline{{\omega }_{\mathbf{K}}^2}=\left(\begin{array}{cc}{\Phi }^{\mathrm{AA}}& {\Phi }^{\mathrm{AB}}\\ {\Phi }^{\mathrm{BA}}& {\Phi }^{\mathrm{BB}}\end{array}\right){\bar{\omega }}_{\mathbf{K}}^2$$

$$\bar{\omega}_{\mathbf{K}}^{2}=\frac{N_{c}}{N}\sum_{\tilde{k}}\sin^{2}\bigg{(}\frac{(K_{x}+\tilde{k}_{x})a}{2}\bigg{)}\\ +\sin^{2}\bigg{(}\frac{(K_{y}+\tilde{k}_{y})a}{2}\bigg{)}+\sin^{2}\bigg{(}\frac{(K_{z}+\tilde{k}_{z})a}{2}\bigg{)}\\$$

$${\Delta }^{\prime }(l,l^{\prime })=\{\begin{array}{ll}{\sum }_{\mathbf{K}}[{\Delta }^{\mathrm{AA}}(\mathbf{K},\omega )]e^{i\mathbf{K}.(R_l-R_{l^{\prime }})},& \text{if l∈A, l′∈A }\\ {\sum }_{\mathbf{K}}[{\Delta }^{\mathrm{AB}}(\mathbf{K},\omega )]e^{i\mathbf{K}.(R_l-R_{l^{\prime }})},& \text{if l∈A, l′∈B }.\\ {\sum }_{\mathbf{K}}[{\Delta }^{\mathrm{BA}}(\mathbf{K},\omega )]e^{i\mathbf{K}.(R_l-R_{l^{\prime }})},& \text{if l∈B, l′∈A }.\\ {\sum }_{\mathbf{K}}[{\Delta }^{\mathrm{BB}}(\mathbf{K},\omega )]e^{i\mathbf{K}.(R_l-R_{l^{\prime }})},& \text{if l∈B, l′∈B }.\end{array}$$

$$D^c(l,l^{\prime },\omega ,V)=$$

$$\sqrt{1-(\hat{V}{)}_l}({\omega }^{2}I-\underline{{\Phi }^{\prime }}-\underline{{\Delta }^{\prime }}-\underline{\hat{V}}{)}_{l{l}^{\prime }}^{-1}\sqrt{1-(\hat{V}{)}_{l^{\prime }}}.$$

$$\displaystyle\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm AA}=\bigg{\langle}$$

$$\displaystyle\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm AB}=\bigg{\langle}$$

$$\displaystyle\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm BA}=\bigg{\langle}$$

$$\displaystyle\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm BB}=\bigg{\langle}$$

$$\underline{{D^{c}_{DCA}}}=\begin{pmatrix}\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm AA}&\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm AB}\\ \\ \big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm BA}&\big{(}D^{c}_{\scriptscriptstyle{\rm DCA}}\big{)}_{\rm BB}\\ \end{pmatrix}$$

Full text

1 Introduction

Most forms of disorder in a crystal structure have two essential consequences, namely a randomness in mass, and a concomitant change in the bond strengths. The phonon spectrum is, naturally, affected strongly by the presence of mass and bond disorder, and can exhibit Anderson localization (AL) Anderson (1958) depending on the nature and strength of disorder, dimensionality and other factors. The possibility of phonon localization due to disorder has evoked great interest over the past several decades, and several theoretical and experimental investigations have been carried out. In recent years, there has been a resurgence of interest in the field due to several direct experimental observations of phonon localization Magdău and Ackland (2017); Howie et al. (2014). The role of AL in the formation of polar nanoregions in ferroelectrics has also generated a lively debate Manley et al. (2014). Theoretical investigations of spectral dynamics in mass and force constant (FC)111We use the term force-constant disorder to imply that the force constants have been randomly chosen and are fixed i.e. no attempt has been made to re-compute force-constants based on structural reorganization. disordered systems received a big impetus with the development of mean-field based approaches Kaplan and Mostoller (1974); Kaplan and Gray (1976); Gray and Kaplan (1981); Grunewald (1976). However, single-site theories are, by construction, incapable of incorporating the full non-local nature of force constants. Additionally, some of the Green’s function based attempts failed to maintain Herglotz analytic properties Mills and Ratanavararaksa (1978); Mookerjee (1973), which are essential to produce physically acceptable results. Various other extensions of single-site theoriesYussouff and Mookerjee (1984); Mookerjee and Singh (1985, 1988); Ghosh et al. (2002) have been attempted. Nevertheless, these perturbative methods are plagued by uncontrollable approximations, and hence are unable to treat AL properly.

As a non-perturbative route to understand AL, exact diagonalization (ED) Payton and Visscher (1968) is the most heavily employed method. Although the method does not suffer from approximations, and yields the disorder averaged spectrum, it does have quite a few disadvantages. The state space increases exponentially with system size implying a severe difficulty in simulating large system sizes. The restriction on system sizes, in turn, leads to difficulties in obtaining information of ALPayton and Visscher (1968). The consideration of three independent force constants namely $\phi_{AA},\phi_{AB}$ and $\phi_{BB}$ is necessary for a minimal description of FC disorder in a binary alloy system. Such a consideration compounds the computational expense involved in ED calculations. Another important exact method for studying AL of phonons is the transfer matrix method (TMM)Pinski et al. (2012). Disorder effects in masses and force constants are intertwined with each other in most of the binary alloys, but the state-of-the-art TMM calculations Pinski et al. (2012) have, thus far, treated the ionic masses and force constants as uncorrelated variables, which is quite unrealistic. In addition, such an approximation can lead to a violation of sum rules. Thus, despite extensive attempts, a satisfactory, reliable method for studying phonon localization in a strong mass and FC disordered binary alloy is still lacking, which calls for further theoretical development.

In a recent workMondal et al. (2017), we described the development and application of a typical medium dynamical cluster approximation (TMDCA) for investigating the effect of mass disorder on the AL of phonons. In this work, we incorporate the effects of FC disorder into the existing framework, thus taking a step closer to realistic disorder. The present study has two objectives, namely (i) the development of a formalism for mass and FC disordered systems, that is non-perturbative, systematically convergent, causal, computationally feasible, and is well-benchmarked; and (ii) the application of this formalism to address several open questions (see below). Our first objective has two parts: to develop (a) the DCA to obtain phonon spectra and (b) the TMDCA to investigate phonon localization, in a mass and FC disordered binary alloy. The DCA yields the average density of states (ADOS), which is experimentally observable and is crucial for a basic understanding of disordered lattice vibrations, at a dramatically less computational cost compared to other methods like ED. Concomitantly, the most salient feature of the TMDCA compared to other theories of localization is its ability in predicting localization based on a single-particle order parameter, namely the typical density of states (TDOS). The main development in this work is a cluster adaptation of the Blackman, Esterling and Berk (BEB) formalism Blackman et al. (1971), which was originally proposed for bond-disordered electronic systems. In this study, we adopt a scalar binary alloy model, that consists of a single branch and a single basis atom within the harmonic approximation.

In order to assess the validity of the method, we carry out quantitative benchmarks, and find excellent agreement. Subsequently, through an application of the method, we attempt to address the following questions/issues related to the effects of pure FC disorder and the interplay of mass and FC disorder: (1) Structurally disordered glasses are known to exhibit a low frequency anomaly, known as the Boson peak, in the density of states. Such an anomaly has also been observed in disordered lattice models, though the origin of the Boson peak in the two families of systems might be quite different. It has been argued that a purely FC disordered system can also exhibit such behaviour. We carry out a comprehensive analysis of the deviations from low frequency Debye behaviour, and ask- What are the reasons and conditions for the emergence of a Boson peak in a binary mass and spring disordered system? (2) In a related context, are the modes associated with the Boson peak localized or delocalized? (3) As shown in our previous workMondal et al. (2017), lighter isotopic impurities lead to strongly localized, short-wavelength phonons in impurity bands. In this work, we ask whether the effect of FC disorder is to reinforce or negate the localization induced by mass disorder? (4) And finally, as relevant for improving the figure of merit in thermoelectrics, we explore the efficacy of vacancies in realizing strong phonon localization over a broad range of frequencies.

The manuscript is organized as follows. In Section 2, we present our model containing both, the mass and FC disorder, and give a detailed description of the DCA and TMDCA formalism as well as their numerical implementation. In Section 3, we present our results and discussion. We summarize our analysis with some future perspectives in Section 4.

2 Method

The Hamiltonian for lattice vibrations involving a single basis atom is:

[TABLE]

where $u_{\alpha}(l)$ is the displacement of an ion in the $\alpha^{\rm th}$ direction having mass $M(l)$ in the $l^{th}$ unit cell coupled by the force constant $\Phi_{\alpha\beta}(l,l^{\prime})$ tensor to $u_{\beta}(l^{\prime})$. We note that the mass $M(l)$ in Eq.(1) can vary randomly from site to site. Since we are considering the Hamiltonian in Eq.(1) for a binary alloy, the site $l$ can be occupied by either an A-type atom or a B-type atom, i.e $M(l)\in\{M_{A},M_{B}\}$ with certain probabilities depending on the relative concentrations of the A or B-type atoms. To this end, it is convenient to introduce occupation indices $(x,y)$ for host (A-type atoms) and guest (B-type atoms) as (following Blackman, Ester and Berk (BEB) Blackman et al. (1971))

[TABLE]

These occupation indices must obey the following properties:

[TABLE]

Note that double occupancy of a a given site is prohibited in this formalism. With this assumption, we are ready to express the randomness in the masses as

[TABLE]

We incorporate such randomness in our formalism by defining a local disorder potential matrix ${\hat{V}}$, as

[TABLE]

where $M_{0}$ is a reference (host) mass, and as a convention, we have chosen A to be the host, hence $M_{0}=M_{A}$.

The corresponding probability distribution for binary disorder reads as

[TABLE]

where $c_{\rm A}$ and $c_{\rm B}=1-c_{\rm A}$ are the concentrations of A and B type of atoms respectively.

Mass disorder can be isotopic or non-isotopic. In general, mass disorder will be accompanied by a corresponding randomness in the force constants as

[TABLE]

Note that $\Phi(l,l^{\prime})$ can be decomposed into diagonal $\Phi(l,l)$ and off-diagonal parts

[TABLE]

where $\Phi_{D}$ and $\Phi_{nn}$ are the diagonal, and the off-diagonal component of the tensor, respectively, and $\vec{\delta}$ is defined as a vector from a site to its nearest neighbors. Force constant tensor must obey a sum rule, namely $\sum_{l^{\prime}}\Phi^{\alpha\beta}(l,l^{\prime})=0$. For satisfying the sum rule, the formalism must incorporate multi-site correlations, because the spring constant tensor is off-diagonal in nature. We are satisfying this property by systematically by increasing $N_{c}$.

Now, we are in a position to apply the equation of motion method to the Hamiltonian (Eq.(1)) and obtain the Dyson equation as

[TABLE]

Next, we premultiply and postmultiply the above Eq.(9) by $x_{l}$ and $y_{l}$, which would generate the four possible configurations of the binary alloy. Combining this with Eqs.(3), (4) and (7) yields four self-consistent equations for the Green’s functions as given below:

[TABLE]

where the four configuration dependent Green’s function are defined as

[TABLE]

and the bare locators, $d_{\rm A}$ and $d_{\rm B}$ are given by

[TABLE]

The four self-consistent equations, Eqs.(10)-Eqs.(13) may be combined in a convenient $2\times 2$ matrix form as

[TABLE]

[TABLE]

And finally, even this matrix equation can be compactified to get

[TABLE]

Here $\underline{D}$ and $\underline{d}$ are matrices of size $2N\times 2N$, where $N$ is the system size.

Thus, we have obtained an equation which has a structure similar to the one obtained in the BEB formalism Blackman et al. (1971) for the electronic problem. It is interesting to note that there is no randomness associated with $\Phi$ matrices. All the randomness is absorbed in the $d$ matrices, and the origin of this randomness lies in the mass term. The $\Phi$ matrices will take the values depending on the random values associated with the mass term. Hence, diagonal mass disorder and off-diagonal spring disorder are dependent on each other. The other point to note is that we can consider three different spring constants, which has been a computational limitation for some theoretical approaches Pinski and Roemer (2011); Beltukov et al. (2013). In order to solve these equations, we have adopted the dynamical cluster approximation. The formalism is similar to the one presented in our previous work on mass disorder, but there are certain steps that are unique to the spring disorder case. In the next section, we provide the details of the formalism.

2.1 Dynamical Cluster Approximation (DCA)

The advent of dynamical mean field theory (DMFT) led to a sort of revolution in the understanding of quantum many body lattice systems. However, since DMFT ignores non-local dynamical correlations, several phenomena such as d-wave superconductivity, Anderson localization, and low dimensional physics are out of scope of this framework. Hence, quantum cluster approaches, that go beyond DMFT, have assumed great importance. One such approach, that is based on momentum space clusters, is the dynamical cluster approximation (DCA).

The DCA may be viewed as an approximation to the wave-vector sums that occur in Feynman-Dyson Perturbation Theory. Here, the first Brillouin zone containing $N$ wavenumbers ${\mathbf{k}}$ is broken into $N_{c}$ non-overlapping coarse graining cells. We then approximate the integrals associated with each diagram by its sum of average/coarse-grained estimates of the integrand within each cell. There is considerable freedom in how this is done. For example, if the integrand is composed of the product of two functions of the integration variable, do we take the product of the two averages, or the average of the products to define the approximate value in the cell? Since these two approximations have the same error provided that the number of such cells is large, we can use this freedom to simplify the approximation. To do this, we define the many-to-few mapping $M({\mathbf{k}})={\mathbf{K}}$, where ${\mathbf{K}}$ labels the cells including ${\mathbf{k}}$, so that ${\mathbf{k}}={\mathbf{K}}+{\mathbf{\tilde{k}}}$, where ${\mathbf{\tilde{k}}}$ labels the wave-numbers within each cell. The corresponding transformation of the Lie algebra is ${\bar{c}}_{\mathbf{k}}=\sum_{\mathbf{\tilde{k}}}c_{{\mathbf{K}}+{\mathbf{\tilde{k}}}}\,.$ It is easy to see that this transformation preserves the Lie algebra. So it maps bosons onto bosons and fermions onto fermions. The mapping is not canonical though since information is lost in the process. Nevertheless, this mapping ensures that FDPT may be used to analyze lattice plus the quantum impurity problem. Under this transformation, all points within each cell are considered to be equivalent, and are mapped to a single point ${\mathbf{K}}$. So, each Green’s function within the cell may be replaced by its average value within the cell. Equivalently, each $G({\mathbf{k}})$ in a Feynman graph may be replaced by its coarse grained analog. More significantly, each sum over ${\mathbf{k}}$ is replaced by a sum over ${\mathbf{K}}$ thereby dramatically reducing the complexity of the problem of order $N$ to order $N_{c}$. The associated FDPT is the same as a small self-consistently embedded periodic cluster problem. Once the cluster problem is solved, we calculate the corresponding irreducible self energy and vertex functions. We use them in the Dyson and Bethe-Salpeter equation to calculate the single-particle spectra and the two-particle susceptibilities.

For example, the DMFT framework may be represented as a mapping of the entire first Brillouin zone to just one momentum at the centre for the zone. The main simplification in the DMFT framework is the absence of momentum conservation at the vertices in the Feynman diagrams, thus leading to a local self-energy. The DCA targets this lacuna of DMFT and replaces the Dirac delta function that represents true momentum conservation at the vertices by a Laue function that conserves momentum, but only for the cluster momenta. This brings back the momentum dependence in the Green’s functions and self-energy, lost at the DMFT level. Thus, as the number of clusters increases, the Brillouin zone is sampled more densely, and hence the thermodynamic limit is systematically approached. We refer the reader to review articles Maier et al. (2005); Terletska et al. (2018) for details of the DCA, and its applications to a variety of problems. The DCA algorithm that we have implemented is derived using the formalism described in the previous subsection and is described below:

1. We start with an initial guess of the hybridization function as

[TABLE]

This guess may be obtained through a coarse graining of the non-disordered Green’s function, or from a previously converged calculation.

2. As a first step for solving the cluster problem, we generate random configurations of the disorder potential V. The disorder potentials $V_{A}$ and $V_{B}$ are assigned depending on whether the site is occupied by A or B type of atom. We generate some random number and if it is less than a given impurity concentration $c_{A}$, we assign a given site as A-type, else it is assigned as B-type.

3. We define $\Phi^{\prime}$ as configuration dependent force-constants which can be obtained by configuration dependent Fourier transform as shown below

[TABLE]

where the 2X2 dispersion matrix is given by

[TABLE]

and the coarse-grained dispersion is given by

[TABLE]

The real space hybridization function $\Delta^{\prime}$ is obtained from the configuration dependent Fourier transform as below:

[TABLE]

After constructing $\Phi^{\prime}$, $\Delta^{\prime}$ and $V$, we compute the corresponding cluster Green function through the mass-weighted Dyson’s equation Mondal et al. (2017) as

[TABLE]

4. The next step is disorder averaging over disorder configurations denoted by $\langle(...)\rangle$. These disorder averaged Green’s functions correspond to a translationally invariant system, and are denoted by the DCA subscript:

[TABLE]

Next, we construct a matrix of the cluster Green function by re-expanding the Green function to a $2N_{c}\times 2N_{c}$ matrix. It can be represented as

[TABLE]

5. As mentioned above, after disorder averaging, the translation symmetry is restored and we can perform Fourier transform for each component to get disorder averaged $\mathbf{K}$ dependent cluster Green function as

[TABLE]

6. Once the cluster problem is solved, we calculate the coarse-grained lattice Green function as

[TABLE]

which in explicit matrix form is given by,

[TABLE]

7. The DCA self consistency condition requires that the disorder averaged cluster Green function equal the coarse-grained lattice Green’s function

[TABLE]

8. The self consistency condition is used for updating the hybridization function for each component

[TABLE]

In the above equations, the $\xi$ is a mixing parameter that determines the fraction of the updated hybridization that should be mixed with the existing one, thus ensuring smooth convergence of the DCA iterations.

As is well known, the Anderson localization of phonons requires us to go beyond DCA. The arithmetic averaging procedure needs to be modified, and a typical averaging ansatz needs to be evolved. Such an ansatz has been worked out in the electronic case, and has been benchmarked against known results. Terletska et al. (2014). We have adopted the same ansatz in the phonon case, and have found that it yields the same level of benchmarks as the electronic case. The formalism that employs this typical averaging ansatz is called the typical medium DCA or the TMDCA, and is detailed in the next section.

2.2 Typical Medium Dynamical Cluster Approximation (TMDCA)

As mentioned in our previous discussionMondal et al. (2017), the effective medium is constructed via algebraic averaging in the DCA, while the TMDCA utilizes geometric averaging to construct the effective medium. We employ the same ansatz for evaluating the typical density of states as in the electronic case Terletska et al. (2014), as:

[TABLE]

where $D^{c}_{ii}$ is defined as

[TABLE]

Next, we calculate the cluster-averaged typical Green function which is also a $2\times 2$ matrix

[TABLE]

We compute each component of the cluster-averaged typical Green’s function from the corresponding component of the typical density of states (26) using the Hilbert transform,

[TABLE]

Once the disorder-averaged cluster Green function is calculated using (30), the self-consistency follows the same steps as in the DCA presented in previous section. The coarse-grained lattice Green’s function is then calculated using (20), which is utilized to update the hybridization function in (25).

Using DCA and TMDCA, we can calculate the arithmetically averaged density of states (ADOS) and typical density of states (TDOS) respectively as follows:

[TABLE]

where $\sigma,\sigma^{\prime}=A/B$. The formalism described above has been implemented, and we present results and discussion in the following section.

3 Results and discussions

We begin our discussion with a validation of the method. To that end we compute the ADOS via exact diagonalization (ED) of a large number of phonon models of large sized disordered supercells. Within the supercells the impurities are randomly distributed and beyond the supercell boundaries the impurity distributions periodically repeat. Specifically for each set of model parameters we derive the force constant matrices of 100 supercells each with 60 impurities and roughly 400 sites on average. The dynamical matrix of each supercell is then evaluated and diagonalized on a $10\times 10\times 10$ supercell momentum space grid. We note that not only the impurity distributions but also the shapes of the supercells are randomized under the following constraints. The supercell volumes vary within 375 and 430 sites and the angles between the vectors that span the supercells vary between 75 and 105 degrees.

The left panels of figure 1 show DCA results (red solid lines) and ED results (black solid lines) for a binary alloy ($c=0.15;V=0.67$) with three different spring constant combinations. A good agreement is seen over all scales thus validating the formalism. We also note that the DCA can access detailed information of the ADOS with relatively small cluster sizes, i.e $N_{c}=64$ which shows that the DCA is dramatically less expensive compared to ED while being numerically exact. To check the sensitivity of the DCA results on the choice of cluster size $N_{c}$, we present the ADOS for different $N_{c}$ in the right panel of Fig.1. We find that the ADOS for $N_{c}=64$ and $N_{c}=125$ are almost identical which implies a rapid convergence of our calculations with respect to increasing $N_{c}$. In contrast to the results for $N_{c}=64$ and $N_{c}=125$, the single-site ($N_{c}=1$) calculations (solid blue lines) are unable to capture non-local fluctuations and disagree significantly with the converged spectra.

The result presented above, namely a comparison of DCA with ED, lends strong credence to results from DCA. Hence, we employ DCA, and subsequently TMDCA to investigate the interplay of spring and mass disorder on phonon spectra and on AL of phonons. We begin with an investigation of the effect of pure spring disorder on phonon spectra.

The upper panel in figure 2 shows the average DOS for pure spring disorder (i.e. $V\rightarrow 0$) with spring constant values $\phi_{AA}=1.0,\;\phi_{BB}=0.1$ and $\phi_{AB}=0.3$ for various impurity concentrations ($c_{B}$) ranging from 0.95 to 0.05. The parameters have been chosen to mimic the values obtained in a recent experiment Gratale et al. (2013) on crystals of binary hard-soft microgel particles with three distinct interparticle potentials. The spring constant values imply that A are hard particles, while B are soft. Hence, $c_{B}=0.95$ corresponds to B-particle concentration of 95% which implies hard sphere concentration of 5%. As expected, the spectrum for a higher concentration of hard particles (stiffer springs, $c_{B}=0.05$, $c_{A}=0.95,\phi_{AA}/\phi_{BB}=10$) has almost entire spectral weight at higher frequencies, and as $c_{B}$ varies from $0.05$ to $0.95$, spectral weight is transferred to lower frequencies.

The DOS at 20% to 40% soft particles ($c_{B}=0.2-0.4$) shows a clear excess density of states around a frequency, that occurs far below the van Hove singularities of the pure hard particle system. Such behavior is strongly reminiscent of disordered systems, where the origin of such an excess of DOS, termed a Boson peak, has generated a lot of debate. We briefly review a few theoretical and experimental results relevant to this issue, and place our results in perspective.

It has been shown theoretically Schirmacher et al. (1998), that a strongly disordered three-dimensional system of coupled harmonic oscillators with a continuous force constant distribution exhibits an excess low-frequency DOS (boson peak) as a generic feature. Specifically, if the system is proximal to the borderline of stability, a low-frequency peak (i.e the Boson peak) appears in the quantity $g(\omega)/\omega^{2}$ as a precursor of the instability. Our results have been obtained for a binary alloy with three values of spring constants, and we see that a Boson peak appears in a regime of lower soft particle concentration.

Experimental measurements Kaya et al. (2010) of normal modes and the DOS in a disordered colloidal crystal showed Debye-like behavior at low energies and an excess of modes, or Boson peak, at higher energies. The normal modes took the form of plane waves, that hybridized with localized short wavelength features in the Debye regime but lost both longitudinal and transverse plane-wave character at a common energy near the Boson peak. More recently, experimentsGratale et al. (2013) on deformable microgel colloidal particles with random stiffness appear to contradict the theoretical results of RefSchirmacher et al. (1998). The authors create crystals of binary hard-soft microgel particles with three distinct interparticle potentials distributed randomly on a two-dimensional triangular lattice. The nearest-neighbor bonds are either very stiff ($\phi_{AA}$), very soft ($\phi_{BB}$), or of intermediate stiffness ($\phi_{AB}$). Subsequently, they obtain, experimentally, the phonon modes in crystals with bond strength disorder as a function of increasing dopant concentration. The interesting feature of the microgel crystal is that although the bonds are randomly distributed, the masses are nearly identical, hence the disorder is purely off-diagonal. The experimental resultsGratale et al. (2013) show the absence of a boson peak, although an excess in the density of states as compared to conventional Debye behaviour was observed.

In the lower panel of figure 2, we present integrated DOS as a function of frequency. The results indicate conventional Debye behaviour ($\sim\omega^{d}$) at lowest frequencies followed by a deviation, and finally a convergence to the normalization value of one at the highest frequencies. In the experiment, the hard particle concentration has been varied from 0 to about 21%. A comparison to figure 3 of Yodh et. alGratale et al. (2013) shows that our results concur well with the experiments. An absence of Boson peak, as concluded in the experiments, is natural since a clear Boson peak occurs only in the opposite limit of lower soft particle concentration. To summarize, within the framework of DCA, we find (see figure 2), for a binary alloy, that a transfer of spectral weight to lower frequencies results in the Boson peak, which emerges as a crossover feature between a pure host system and a pure guest system.

A question that has been much debated in the literature is about the nature of states within the boson peak - Are they localized or delocalized? This question can be effectively answered through the evaluation of the typical DOS, since the typical spectral weight is a measure of the proximity to Anderson localization. A subtle issue about the interpretation of the typical density of states must be mentioned here. A non-zero typical DOS signifies the presence of extended states. According to Mott, a degeneracy of localized and extended states should lead to their hybridization, and hence an eventual delocalization of the localized states. The average DOS and the typical DOS, being different at a given energy, is thus immaterial regarding the identification of the states being extended or localized. If the typical DOS is non-zero, the states at that energy should be interpreted as being extended. Concomitantly, a large difference between the average and typical DOS does indicate a proximity to the Anderson localization transition (ALT).

In figure 3, we show, for the same parameters as figure 2, a series of average (black solid lines) and the corresponding typical spectra (red shaded part).

The boson peak is seen to have numerically negligible typical spectral weight (TSW) until $c_{B}\sim 0.4$, beyond which the low frequency peak acquires finite and significant TSW. Thus, the states in the BP exhibit a kind of crossover from being almost localized (proximal to ALT) to being delocalized (relatively smaller difference between average and typical DOS) with increasing $c_{B}$. The overall typical spectral weight, shown in the bottom left panel, is non-monotonic, and shows a minimum at $c_{B}\sim 0.5$, showing that the overall system is closest to the ALT when the ratio of the concentrations of the two species is roughly equal to one.

In order to understand the interplay of mass and spring disorder, we consider two protocols. In the first, we keep the mass ratio parameter, $V=0.67$ and the impurity concentration $c=0.5$ as constants, and increase the spring disorder systematically by varying $\phi_{AA}/\phi_{BB}$ (with $\phi_{BB}=1$, and $\phi_{AB}=\left(\phi_{AA}+\phi_{BB}\right)/2$) from $0.2$ to $2.0$, representing a change of host spring constants from very soft to very stiff. The resulting spectra (ADOS as solid black lines and TDOS as red shaded part) are shown in figure 4, while the inset shows the integrated typical spectral weight (TSW, solid blue circles) as a function of $\phi_{AA}$. For soft A-springs, the characteristic frequencies of the system must be lower than a pure B-type system, and as the $\phi_{AA}$ is increased, spectral weight in the second, high frequency peak increases, as also the bandwidth of the system. So, nominally, it appears that the system is getting delocalized, as the host springs are made stiffer for a fixed mass disorder. However, the inset shows a decrease in TSW with increasing stiffness of $\phi_{AA}$, which implies that the order parameter for AL is decreasing, and hence the system is moving closer to localization. If we focus on a fixed frequency, say $\omega^{2}=5.0$, then we see that for $\phi_{AA}=0.2$, the ADOS and TDOS are zero, while for $\phi_{AA}=2.0$, both the average and typical DOS are non-zero, suggesting the interpretation that spring disorder is of a delocalizing nature and counters the localization produced by mass disorder alone. However, the order parameter for AL, namely the TSW, decreases sharply. Thus the interplay of mass and spring disorder is quite subtle, and an interpretation of the results need to be done carefully. It must be emphasized here, that the subtlety of this interplay has been uncovered through the application of TMDCA, which is able to produce, simultaneously, the typical and the average DOS.

The second protocol is to vary the mass ratio parameter, $V=1-M_{imp}/M_{host}$, keeping the relative concentrations as well the spring constants fixed. The ADOS and TDOS are shown in figure 5 for $V=0.01$ to $0.9$, implying a systematic decrease in the B-site ionic mass. Again, lighter impurities imply transfer of spectral weight to higher frequencies, and the B-site band appears as a separate feature, which blue shifts significantly with increasing $V$. In parallel with the results of the first protocol, this result lends itself to an interpretation of delocalization of modes at higher frequencies, but the vanishing of the typical density of states shows that the high frequency modes for $V\rightarrow 1$ are almost localized. The inset shows that the TSW (solid blue circles) decreases sharply with increasing $V$, and this also implies that increasing mass disorder in the presence of fixed spring disorder pushes the system closer to the AL transition.

The insight we gain from the study of the interplay of mass and spring disorder is that an inference of localization/delocalization of specific modes cannot be made on the basis of ADOS alone, and the TDOS must be concomitantly examined.

Vacancies, even at low concentrations, can lead to strong localization of phonons. Within the present theoretical framework, we model vacancies as weakly bonded sites with vanishing mass. So, we choose $M_{imp}=M_{host}/20$, which is equivalent to $V=0.95$, and spring constants as $\phi_{AA}=1,\phi_{BB}=0.15,\phi_{AB}=0.15$. For these parameters, in figure 6, we show the average DOS (upper panel) and typical DOS (lower panel) for three different guest concentrations, namely $c=0.1$(solid black), $c=0.2$(dashed red), and $c=0.3$(dotted blue). The upper panel shows that the average DOS hardly changes with increasing concentration, while the corresponding typical DOS (lower panel) undergoes significant changes. The inset in the lower panel shows the rapid decrease of the integrated typical spectral weight (blue solid circles) with increasing concentration, $c$. Modeling real vacancies is quite challenging, but the present analysis with a very crude model for vacancies is already indicative of their strong localization effects. The figure of merit for thermoelectrics is inversely proportional to the thermal conductivity, and directly proportional to electrical conductivity. So, in order to maximize the figure of merit, the vacancy concentration, $c$ should be tuned to an optimal value such that it is less than, but not too close to the percolation limit, implying that the electrical conductivity is not too significantly affected, but the thermal conductivity due to phonons gets drastically reduced due to the strong localization of acoustic phonons in a large part of the spectrum.

Finally, we attempt a qualitative comparison with experiments. It has been argued for a NixPt1-x alloy, that $x=0.65$ constitutes weak force-constant disorder, while $x=0.5$ constitutes strong mass and force-constant disorder. Ghosh et. al.Ghosh et al. (2002) demonstrate that for $x=0.5$, a CPA-level consideration of inter-atomic force-constants leads to a split-band spectrum. The authors show that a proper treatment of force-constant disorder using itinerant CPA leads to a closure of the gap. Our calculations are in full qualitative agreement with these conclusions as argued below. Figure 7 shows ADOS for a binary alloy with Mimp=Mhost/3 as appropriate for Ni impurities in Pt host. The left panel is for $N_{c}=1$, equivalent to a CPA calculation, while the right panel is for $N_{c}=64$, which is equivalent to the thermodynamic limit. The impurity concentration used is $c=0.5$, which implies strong mass disorder; and the spectra corresponding to three distinct force-constant combinations are shown. The black solid line corresponds to pure mass-disorder, which at the CPA-level shows a split-band (left panel), while at the DCA-level (right panel), the spectrum has a two peak structure with a soft-gap between the peaks. The red and green lines correspond to weak and strong force-constant disorder respectively. Again, the CPA results are hardly affected by an increase in disorder, while the DCA results for $N_{c}=64$ show that increasing force-constant disorder leads to significant spectral weight transfer, especially a filling-up of the soft-gap. Since the force constant combination represented by the green line most closely corresponds to the Ni0.5Pt0.5 alloy, we conclude that our results agree qualitatively with the ICPA results as well as experimental neutron scattering data for Ni0.5Pt0.5 Tsunoda et al. (1979).

4 Conclusions

The incorporation of the BEB formalism for off-diagonal disorder into the TMDCA yields a reliable, and computationally feasible approach for investigating binary mass and spring-disordered alloys. Such a conclusion is supported by the benchmarking studies discussed in the initial part of the results section. For a fixed mass ratio, and fixed spring constants ($\phi_{AA},\phi_{BB}$ and $\phi_{AB}$, increasing the soft particle concentration leads to an excess density of states below the first van Hove singularity of the host hard particle system. In the present context, it may be identified as the Boson peak, commonly observed in structurally disordered glasses as well as disordered lattice systems, albeit with different origins. We conclude that the origin of the Boson peak in the disordered binary alloy system is a transfer of spectral weight from the guest to the host system, with a necessary condition being the presence of off-diagonal disorder. We emphasize that with pure mass disorder, even though spectral weight transfer does occur, such a BP does not emerge. Additionally, we find that at very low soft particle concentrations, the states in the BP are completely localized, but some of the states crossover to being extended as their concentration is increased. The BP eventually ceases to be an anomaly, as the soft particle system becomes the host, and the hard particles assume the role of impurities. The interplay of mass and spring disorder is found to be quite subtle. The overall typical spectral weight decreases upon increasing either of the types of disorder, which indicates that there is a co-operative interplay. However, an added clause is that the interpretation of a co-operative or competitive interplay is also frequency selective, since different parts of the spectrum can transform from being localized to delocalized or vice-versa, depending on the protocol. With a crude modeling of vacancies, we suggest that tuning the concentration of vacancies to an optimal level, which is below but not close to the percolation limit, should be an optimal route to maximizing the figure of merit of thermoelectric materials. All the above results finally culminate in an attempt to understand experiments on NixPt1-x alloys, where, in agreement with ICPA results, we show that $x=0.65$ constitutes weak force-constant disorder, while $x=0.5$ represents a system with strong mass and force-constant disorder. Of course, such a conclusion is qualitative at best, because the real system has a non-trivial structure with multiple branches and spring-constant disorder, while the present implementation is restricted to a scalar approximation and force-constant disorder. The extension of the present framework to incorporate multiple branches, which will allow us to treat disordered phonons in real materials, is in progress.

Acknowledgements.

A portion of this research (T.B.) was conducted at the Center for Nanophase Materials Sciences, which is a Department of Energy (DOE) Office of Science User Facility. This material is also based upon work supported by the National Science Foundation under the NSF EPSCoR Cooperative Agreement No. EPS-1003897 with additional support from the Louisiana Board of Regents (M.J., W.R.M.). N.S.V and W.R.M acknowledge funding from JNCASR, India. M.J acknowledges support from DOE grant DE-SC0017861.