Structural origins of electronic conduction in amorphous copper-doped alumina

TL;DR

This study uses ab initio modeling to explore how copper clustering and structural features in amorphous alumina influence electronic conduction, revealing the role of Cu concentration and connectivity.

Contribution

It provides new insights into the structural origins of conduction in amorphous copper-doped alumina through detailed modeling and analysis of Cu clustering, electronic properties, and conduction pathways.

Findings

Cu atoms tend to cluster in alumina.

Higher Cu concentrations lead to metallization by filling the band gap.

Conduction pathways depend explicitly on Cu connectivity.

Abstract

We perform an {\it ab initio} modeling of amorphous copper-doped alumina (a-Al$_2$O$_3$:Cu), a prospective memory material based on resistance switching, and study the structural origin of electronic conduction in this material. We generate molecular dynamics based models of a-Al$_2$O$_3$:Cu at various Cu-concentrations and study the structural, electronic and vibrational properties as a function of Cu-concentration. Cu atoms show a strong tendency to cluster in the alumina host, and metallize the system by filling the band gap uniformly for higher Cu-concentrations. We also study thermal fluctuations of the HOMO-LUMO energy splitting and observe the time evolution of the size of the band gap, which can be expected to have an important impact on the conductivity. We perform a numerical computation of conduction pathways, and show its explicit dependence on Cu connectivity in the host. We present an analysis of ion dynamics and structural aspects of localization of classical normal modes in our models.

Tables

Table 1. Table 1: Initial and final densities of a-Al 2 O 3 :Cu models

Cu content	Mol. Formula	${\rho }_{in}$(g/cc)	${\rho }_f$(g/cc)
0%	(Al2O3)1.00Cu0.00	3.175	3.175
10%	(Al2O3)0.90Cu0.10	3.58	3.75
20%	(Al2O3)0.80Cu0.20	3.78	3.99
30%	(Al2O3)0.70Cu0.30	4.53	4.82

Table 2. Table 2: Average coordination numbers around Cu atoms for 10%, 20% and 30% Cu models.

Cu content(%)	Cu-O	Cu-Cu	Cu-Al
10	1.15	5.1	3.0
20	0.68	6.85	2.45
30	0.48	8.27	1.78

Equations

$${\sigma }_{\mathbf{k}}(\omega )=\frac{2\pi e^2}{3m^2\omega \Omega }\sum _{i,j}\sum _{\alpha }[f({ϵ}_{i,\mathbf{k}})-f({ϵ}_{j,\mathbf{k}})]$$

$${\mid ⟨{\psi }_{j,\mathbf{k}}|p^{\alpha }|{\psi }_{i,\mathbf{k}}⟩\mid }^2\delta ({ϵ}_{j,\mathbf{k}}-{ϵ}_{i,\mathbf{k}}-\hslash \omega )$$

$$g_{ij}(\mathbf{k},\omega )=\frac{2\pi e^2}{3m^2\omega \Omega }[f({ϵ}_{i,\mathbf{k}})-f({ϵ}_{j,\mathbf{k}})]\delta ({ϵ}_{j,\mathbf{k}}-{ϵ}_{i,\mathbf{k}}-\hslash \omega ).$$

$$\sigma =\sum _{i,j,\alpha }{g}_{ij}\int d^{3}x\int d^3{x}^{\prime }[{\psi }_j^{\ast }(x)p^{\alpha }{\psi }_i(x)][{\psi }_i^{\ast }(x^{\prime })p^{\alpha }{\psi }_j(x^{\prime })],$$

$$\sigma \approx h^6\sum _{x,x^{\prime }}\sum _{i,j,\alpha }{g}_{ij}{\xi }_{ji}^{\alpha }(x){\xi }_{ij}^{\alpha }(x^{\prime })$$

$$\Gamma (x,x^{\prime })=\sum _{i,j,\alpha }{g}_{ij}{\xi }_{ji}^{\alpha }(x){\xi }_{ij}^{\alpha }(x^{\prime })$$

$${\sigma }_{dc}=\frac{2\pi e^2\hslash \Omega }{m^2}{\mid D_{{ϵ}_f}\mid }^2{N}^2({ϵ}_f)$$

$$⟨r^2(t){⟩}_{\alpha }=\frac{1}{N_a}\sum _i^{N_{\alpha }}⟨{\mid \overrightarrow{r_i}(t)-\overrightarrow{r_i}(0)\mid }^2⟩$$

$$\text{centering}⟨r^2(t)⟩=A+6Dt\text{@add@centering}$$

$$Z(E)$$

$$Z_{\alpha }(E)$$

Full text

I Introduction

Non-volatile memory devices based on resistive switching characteristics have been studied since the late 1960s Simmons and Verderber (1967). In these devices, application of an external bias potential across an electrolyte changes the electrical conductivity of the electrolyte by changing its structure. This process is reversible and can be performed in the time scale of nanoseconds. Three types of resistive random access memory (RRAM) devices have been studied in detail Waser et al. (2009) and these include RRAM based on oxygen vacancies, RRAM based on thermo-chemical effects and RRAM based on the electro-chemical metallization. The later class of devices are also called conducting bridge random access memory or CBRAM. The CBRAM devices are composed of a thin solid electrolyte layer placed between an oxidizable anode ($eg.$ Cu, Ag or TiN) and an inert cathode ($eg.$ W or Pt). The Cu, in its ionic state, is converted into the conducting “filament” by the applied field: the ions are reduced by electrons flowing from the cathode to leave them in their metallic form, although other counter ions (e.g., OH-) may also be involved in this process Tappertzhofen et al. (2013). With the application of a reverse bias, the connectivity of the cluster can be destroyed, and the device is put into a highly electronically resistive state. The details of the mechanism of CBRAMs have been described elsewhere Kozicki and Barnaby (2016); Dietrich et al. (2007). The performance of CBRAM devices has been studied with several materials as the solid electrolyte which include chalcogenides Kund et al. (2005); Lee et al. (2017), insulating metal oxides Chen et al. (2017); Tsuruoka et al. (2010); Gu et al. (2010); Xu et al. (2014); Belmonte et al. (2013, 2015); Pandey et al. (2015) and bilayer materials Tsai et al. (2016); Barci et al. (2016). CBRAM devices have demonstrated excellent performance in terms of operational voltage, read/write speed, endurance and data retention. Among the host materials reviewed for CBRAM devices, alumina (Al2O3) shows particular promise. It has a high dielectric constant, large band gap, and its amorphous phase is highly stable Kittl et al. (2009); Eklund et al. (2009). The experimental results for CBRAM devices based on Cu alloyed with Al2O3 have shown that the cell exhibits highly controlled set and reset operations, fast pulse programming (10 ns) at low voltage (<3 V) and low-current (10 $\mu$A) with 106 cycles per second for the writing speed Belmonte et al. (2013).

In this paper, we use ab initio molecular dynamics (AIMD) to generate atomic models of a-Al2O3:Cu and investigate the microscopic origin of elctronic conduction in this material. The work presented in this paper shows that an increase in local Cu-concentration can result in stable conducting pathways due to the strong tendency of Cu atoms to cluster in the ionic host. This would lead to a highly stable low resistance state (LRS) for high copper concentration, which does indeed seem to be the case for copper-alumina devices Belmonte et al. (2013). We study the electronic properties for these models and are able to crudely estimate the local concentration of Cu above which CBRAM device switch to the LRS. We present the numerical computation of conduction-active parts of the network based on our recent work on computing space projected conductivity (SPC) Prasai et al. (2018), and show that the strong electron-lattice coupling for electron states near the gap leads to interesting and substantial thermally induced conductivity fluctuations on a picosecond time scale.

The rest of the paper is organized as follows. Section II describes computational details used to create the structures and also the details of our method to obtain the SPC. Section III includes results where we discuss structural, electronic and vibrational properties of the models in different subsections. Section IV provides the conclusions.

II Computations

II.1 Model Generation

In this work, we use AIMD to generate four atomic models with the composition of (a-Al2O3)1-nCun with $n=$ 0, 0.1, 0.2 and 0.3. We used a density of 3.175 g/cm3 for a-Al2O3 Gutiérrez et al. (2000); Vashishta et al. (2008). For the Cu-doped models, we referred to the literature Miranda Hernández J.G. and Rocha-Rangel (2006) to make an initial guess, then carried out a zero-pressure relaxation to correct/optimize the result. For each model, we began by taking a cubic supercell of 200 atoms with randomly initialized positions of the atoms.  Plane wave density functional calculations were performed using the VASP package Kresse and Hafner (1993) and projector-augmented wave (PAW) Blöchl (1994); Kresse and Joubert (1999) potentials within the local density approximation (LDA) Perdew and Zunger (1981) using periodic boundary conditions. We used a kinetic energy cutoff of 420 eV and the $\Gamma$-point to sample the Brillouin zone. A time step of 1.5 fs was used and the temperature was controlled by a Nosé-Hoover thermostat throughout.

We performed a melt-quench simulation Drabold (2009) with a starting temperature of 3500 K. After annealing the “hot liquid” for 7.5 ps at 3500 K, we cooled each model to 2600 K at a rate of 0.27 K/fs as discussed in reference Momida et al. (2006a) and then equilibrated for 10 ps. Each model was then quenched to 300 K at the same cooling rate 0.27 K/fs and further equilibrated for another 10 ps. Zero pressure relaxations were used to determine the final densities for Cu-doped models. The final force between the atoms is no more than 0.01 eV/atom. The initial and final densities are provided in table 1.

II.2 Spatial Projection of Electronic Conductivity

In this section, we discuss a method to obtain a space projected electronic conductivity. We discuss the method in detail in Ref. Prasai et al. (2018). We begin by writing the diagonal elements of the conductivity tensor for each k-point k and frequency $\omega$ using the standard Kubo-Greenwood formula KGF Kubo (1957); Greenwood (1958) as:

[TABLE]

In the above equation (1), e and m represent the charge and mass of the electron respectively. $\Omega$ represents the volume of the supercell. We average over diagonal elements of conductivity tensor($\alpha=x,y,z$). $\psi_{i,\bf k}$ is the Kohn-Sham orbital associated with energy $\epsilon_{i,\bf k}$ and $f(\epsilon_{i,\bf k})$ denotes the Fermi-Dirac weight. $p^{\alpha}$ is the momentum operator along each Cartesian direction $\alpha$. Let

[TABLE]

Then suppressing the explicit dependence of $\sigma$ on k and $\omega$, the conductivity can be expressed as:

[TABLE]

a form that reminds of the the current-current correlation function origins of Kubo’s approach. If we define complex valued functions $\xi_{ij}^{\alpha}(x)=\psi_{i}^{*}(x)p^{\alpha}\psi_{j}(x)$ on a real space grid (call them x) with uniform spacing of width h in three dimensions, then we can approximate the integrals as a sum on the grid. Thus, Eq. (2) can be written as:

[TABLE]

In the preceding, the approximation becomes exact as $h\rightarrow 0$. If we define a Hermitian, positive-semidefinite matrix:

[TABLE]

we can spatially decompose the conductivity at each grid point as $\zeta(x)=\mid\sum_{x^{\prime}}\Gamma(x,x^{\prime})\mid$. $\zeta(x)$ contains vital information about the conduction-active parts of the system111A simpler scheme is to just look at the structure of the Kohn-Sham eigenfunctions near the Fermi level to identify the conduction active parts of the network. While this is a sensible first approximation, it entirely neglects the current-current correlations that underlie the derivation of The Kubo formula from linear response theory..

To implement the method, we used VASP and associated Kohn-Sham orbitals $\psi_{i,\bf k}$. We divided the supercell into $36\times 36\times 36\ (dim\Gamma=46656)$ grid points and obtained the wavefunction at each point by using the convenient code of R. M. Feenstra and M. Widom Feenstra and Widom (2012). In computing the $\xi_{ij}^{\alpha}$, we used a centered finite-difference method to compute the gradient of $\psi_{i}$ for each $\alpha$. We used an electronic temperature of T = 1000 K for the Fermi-Dirac distribution. We approximated the $\delta$ function in Eq. (1) by Gaussian distribution of width kT, where k is Boltzmann’s constant.

III Results

III.1 Bonding and topology of the models

As a test of validity of our models, we compute the total radial distribution function, $g(r)$, on a-Al2O3 models and compare with experimentally measured neutron scattering $g(r)$ from Lamparter and Kniep (1997). A plot showing these two $g(r)$ is presented in Fig. 1 and shows that the models capture the structural order upto 6 Å reasonably well. We also compute the structure factor, $S(q)$, on our models at 2600 K and compare it with $S(q)$ measured on l-Al2O3 Landron et al. (2001). The plot shows that these two $S(q)$ show a satisfactory agreement, especially on the positions of peaks at 1.8 Å*-1*, 2.8 Å*-1*, 4.7 Å*-1*. The bottom left plot in Fig. 1 presents the partial $g(r)$ computed on models of a-Al2O3. The peaks at 1.81 Å, 2.78 Å and 3.17 Å correspond to the geometrical bond distances for Al-O, O-O and Al-Al pairs respectively; these results are in agreement with similar earlier works Chagarov and Kummel (2009); Momida et al. (2006b); Sankaran et al. (2012). The bottom right plot in Fig. 1 shows the partial $S(q)$ corresponding to Al-Al, Al-O and O-O pairs computed on a-Al2O3 models. We see that the first peak in the total $S(q)$ occurs at 2.8 Å*-1* due to the partial cancellation arising from Al-O correlations.

For doped models, the computed $g(r)$ are plotted in Fig. 2 and shows that the position of first peak remains largely the same as undoped a-Al2O3 suggesting that Al-O bond remains unaltered. As the concentration of Cu increases, a hump corresponding to Cu-Cu correlation appears and grows at $r\approx$ 2.44 Å. The relative sharpness of Cu-Cu hump, even for the lowest concentration of Cu, provides a hint that Cu atoms are probably clustered. Indeed, a visual inspection of the models, shown here in Fig. 3, clearly shows the strong tendency of Cu-atoms to cluster.

It is significant that Cu strongly tends to cluster. A study by Dawson and Robertson Dawson and Robertson (2016) asserts that the Cu-Cu interactions become more favorable with increasing Cu content. We study the average coordination number around Cu atom at different Cu-concentrations as shown in table 2. We take the first minima in partial $g(r)$ as the cutoff distance to define the coordination number. The increase in Cu-coordination by Cu and the decrease in Cu-coordination by Al and O supports the segregation of Cu from the host and formation of cluster.

III.2 Electronic structure

III.2.1 Density of States and the Localization

Doping by copper in a-Al2O3 is expected to have effects on electronic properties which are of interest for applications of these materials in CBRAM devices. We investigate these effects by examining the density of Kohn-Sham eigenstates (EDOS) and their spatial localization. The localization is gauged by computing the inverse participation ratio (IPR) that is defined as IPR=${\sum_{i}{a_{ni}}^{4}}/{(\sum_{i}{a_{ni}}^{2})^{2}}$ Ziman (1979), where the $a_{ni}$’s are the contribution to eigenfunction $\psi_{n}$ from the $i^{th}$ atomic projected orbital obtained from VASP. Fig. 4 shows the computed EDOS and IPR as a function of Cu-concentration. We find a decrease in HOMO-LUMO gap with increasing Cu-concentration; at Cu-concentration 20% and 30%, The EDOS is continuous across the Fermi level. The states that fill-in the band gap are quite extended as indicated by small values of IPR around the Fermi level in Fig. 4. The mean IPR values around the gap declines monotonically with Cu-concentration.

By projecting the electronic states onto atomic sites, we observe that the states near the Fermi level for the doped models consist of Cu-orbitals. An example of the site projected EDOS, for 20% Cu, is plotted in Fig. 5. It is quite interesting that at 20% and 30% Cu-concentrations, Cu levels almost *uniformly *fill the host a-Al2O3 gap. The Cu does not form an impurity band, as one might naively suppose from experience on heavily-doped semiconductors. We see that models with higher Cu-concentration produce states near Fermi level that yield an essentially metallic form of conduction. This is qualitatively different than the case of Ag in GeSe3 Prasai et al. (2017), wherein the Ag atoms do not cluster and do not introduce states in the optical gap of the host. We observe that electron states in the gap are filled mostly by 3d, 4s and 4p orbitals of Cu.

III.2.2 Charge analysis on Cu atoms

The formation of Cu-cluster in a-Al2O3 matrix leaves the Cu atoms in different charge states depending on the local environment of these Cu atoms with O and/or Al atoms. We performed Bader charge analysis Tang et al. (2009) to calculate net charge on these atoms and an analysis for 20% Cu-doped model is shown in Fig. 6. The charge state of the Cu atoms (shown in color in Fig. 6) can be explained by a simple analysis of the first neighbors around the Cu atoms. Among all the Cu-atoms shown in the figure, only five Cu atoms have exclusively Cu neighbors and are neutral in nature; the rest of the Cu are neighbors with at least one Al or O atoms. When a Cu atom is a neighbor with Al or O atoms, bonding or charge transfer occurs. A Cu atom bonded with O atoms is positively charged, whereas a Cu atom bonded with Al atoms is slightly negatively charged and can be understood in terms of difference in electronegativities of Cu and Al. When a Cu atom is bonded with both O and Al atoms, it is charge neutral. The charge compensation likely to happen in such bonding. The Cu atoms shown in green are therefore almost metallic in nature and are likely to form a conducting channel for the current to flow in the network.

III.2.3 Thermally driven conduction fluctuations

In this section, we discuss relatively dramatic thermally-induced fluctuations in the HOMO-LUMO splitting and consider the electronic conduction mechanisms222Here and elsewhere in this paper, electronic time evolution refers only to variation in Kohn-Sham eigenvalues/states on the Born-Oppenheimer surface – no attempt is made to solve a time-dependent Kohn-Sham equation. We illustrate with one of the conducting models (including 20% Cu) and performed molecular dynamics (MD) at 1000 K for 24 ps. The fluctuation of the frontier HOMO and LUMO levels with time is provided in Fig. 7. $\eta(t)$ is the HOMO-LUMO splitting through the course of the MD. The model reveals a large thermally driven fluctuation in the value of the HOMO-LUMO gap with time.

To physically interpret the connection of the gap with electronic conductivity ($\sigma$), let us write a simple expression for the dc conductivity (T = 0 K) following Mott and Davis Mott and Davis (1979),

[TABLE]

where $D_{\epsilon_{f}}$ is a matrix element of $\nabla_{\alpha}$ between Kohn-Sham states near the Fermi level and $N$($\epsilon_{f}$) is the density of states. So, for dc conduction to occur, there needs to be finite density of states at the Fermi level (to enable electronic transitions, as from Fermi’s Golden Rule) and non-vanishing matrix elements ${\mid D_{\epsilon_{f}}\mid}^{2}$ as in Eq. (1). We expect more available states near the Fermi level for the system with small gap, thus the conductivity $\sigma(t)$ can be very crudely linked to $\eta(t)$ (small $\eta\implies$ large $\sigma$) in the spirit of Landau-Zener tunneling Landau (1932); Zener (1932). We provisionally interpret the small gap (small $\eta$) instantaneous configurations as low resistance states, and the large gap configurations as high resistance states.

It is therefore interesting to visualize the conduction-active parts of the network for these different states. We selected two snapshots (shown by orange arrows in Fig. 7), one representing a small gap (low $\eta$) and the other large gap (high $\eta$) from the simulation and obtained the SPC as described in section II.2. The variation of the HOMO-LUMO gap due to thermal fluctuations has also been studied in Boron-doped $a$-Si at 600 K, where it was observed that with addition of hydrogen to the network, there occurs a thermal modulation of HOMO and LUMO states causing the HOMO and LUMO states to be overlapped at a certain interval of the thermal simulation representing highly conducting configuration Pandey et al. (2014). This computation makes it clear that the DC conductivity is difficult to accurately estimate, since to handle the large electron-phonon coupling for states near the Fermi level, long MD averages at constant temperature would be required (within an adiabatic picture for which one simply averages the Kubo formula over a trajectory.

III.2.4 Space-Projected Conductivity

We investigated SPC by computing $\zeta(x)$ as described in section II.2 in our models. SPC values are evaluated at coarse 3D grid points inside the supercells. A graphical representation of SPC values in 3D grid points overlaid with the atomic configuration is shown in Fig. 8. This figure shows the SPC computed on two models: one with large $\eta$ and the other with and small $\eta$. We include 12% of the highest local contributions to SPC in each plot. The SPC reveals that the conduction path is primarily along interconnected Cu atoms. A few O atoms in the vicinity of Cu atoms also participate in the conduction whereas Al atoms do not show any role in the conduction. We see that the SPC for the large gap snapshot is disconnected so that $\zeta(x)$ appears to be localized in certain region whereas the SPC with small gap forms an interconnected chain for the conduction. For these two particular structures, we observed the local configurations as shown by the enclosed circles of Fig. 8 where the Cu atoms come closer to form short bonds and form a closed network. This shows that the connectivity among Cu atoms determines the conductivity of the system. Besides the structural difference, the type and the number of clusters also affect the HOMO-LUMO gap. It has been shown that an alternation of the HOMO-LUMO gap occurs between even and odd numbered isolated clusters due to electron-pairing effects and particularly large gap for cluster size 2, 8, 18, 20, 34 and 40 which are also called as magic clusters Kabir et al. (2004). At this temperature, the diffusion of Cu atoms may cause the change in the bonding environment of Cu atoms resulting in the variation of the gap with time.

III.3 Ionic motion

As a representative example, the 20% model was annealed at different temperatures 800 K and 1000 K for 15 ps, and the resulting ion dynamics were studied by calculating the mean-squared displacement for each atomic species as:

[TABLE]

where $N_{\alpha}$ represents the number of atoms of species $\alpha$, $r_{i}(t)$ represents the position of atom $i$ at time t, and the $\langle$ $\rangle$ represents an average on the time steps and/or the particles. The connection between mean-squared displacement and the self-diffusion coefficient is given by Einstein’s relation

[TABLE]

where $D$ is the self-diffusion coefficient, $A$ is a constant and $t$ is the simulation time.

Figure 9 shows the mean-squared displacement for the corresponding species. Clearly, Cu atoms are more diffusive than Al and O atoms. On taking the snapshots of the position of atoms (figures not shown here), we find that the Cu atoms do not diffuse into the host matrix but diffuse within the Cu clusters and thus the Cu clusters become stable at these range of temperatures. We then calculated the self-diffusion coefficient for each species using Eq. (7). The diffusion coefficient for Cu at 800 K and 1000 K are obtained to be $9.95\times 10^{-7}cm{{}^{2}}/s$ and $6.248\times 10^{-6}cm{{}^{2}}/s$ respectively. Cu is relatively static in a-Al2O3 compared to chalcogenides Prasai and Drabold (2011).

III.4 Lattice Dynamics

We study the lattice dynamics of these Cu-doped systems by the means of vibrational density of states (VDOS), species projected VDOS and the vibrational IPR. The properties are studied within the harmonic approximation using the first principles method. The dynamical matrix is obtained by displacing each atoms by 0.015 Å along $\pm x$, $\pm y$ and $\pm z$ directions. The diagonalization of the dynamical matrix yields eigenfrequencies and the corresponding eigenmodes. The normalized VDOS and the partial VDOS are expressed as Pasquarello et al. (1998)

[TABLE]

where $\omega_{n}$ are the normalized eigenfrequencies (3N in total). Here, the sum over $i$ is over all the atoms belonging to the species $\alpha$ and $e_{i}^{n}$ corresponds to the displacement vector of atom $i$ with Cartesian components $e_{i\mu}^{n}$ where $\mu$ = x, y and z. We approximate the $\delta$ function by a Gaussian distribution function of width 10 $cm^{-1}$. Among the 3N eigenmodes, we neglect the first three translational modes with frequency very close to zero.

Figure 10 shows the total and partial VDOS for 10% and 30% Cu content. The lower vibrational modes correspond to the Cu atoms. The higher frequency modes are unsurprisingly dominated by O atoms. To study the localization of the vibrational eigenstates, we calculated the vibrational IPR for each species. From Fig. 11, we see that the higher modes corresponding to the O atoms are more localized compared to the lower modes for both concentrations of Cu. The lower eigenstates corresponding to Cu for 10% Cu model are quite localized compared to the 30% Cu model. The vibrational states for aluminum are mostly extended for both models.

IV Conclusion

In this paper, we studied realistic models of a-Al2O3:Cu, and showed that the Cu atoms have a strong propensity to cluster in the ionic a-Al2O3 host. We observed a continuous filling of the optical gap by Cu levels, especially at 20% and 30% models. As the Cu-concentration increases (and Cu-Cu connectivity increases), the Cu levels band to enable metallic conduction. We observed the opening and closing of the HOMO-LUMO gap at an elevated temperature, and projected electronic conductivity into real space and visualized the conduction-active parts of the network. We showed that the connectivity of Cu atoms play a significant role in the electronic conduction. We studied the diffusion of Cu atoms in a-Al2O3 at different temperatures and observed that the Cu atoms do not diffuse easily into the a-Al2O3 in contrast with relatively covalent chalcogenides like GeSe3 Prasai and Drabold (2011). We discussed the harmonic lattice dynamics of the models by calculating vibrational density of states and the vibrational IPR and showed that the lower vibrational modes correspond to Cu atoms and the higher modes correspond to O atoms.

The results presented in this work on a-Al2O3:Cu show an interesting contrast with similar study performed on GeSe3:Ag Subedi et al. (2019). We find that the properties of Cu in the oxide host (in this case, a-Al2O3:Cu) contrast with that of Ag in chalcogenide (in case of Subedi et al. (2019), GeSe3:Ag). The Ag atoms do not form a cluster in the GeSe3 and no uniform filling of the optical gap is observed. In other words, one has to electrochemically work hard to draw Ag atoms together to form a cluster in GeSe3. So, the electronic conduction is likely to occur by hopping process in GeSe3:Ag whereas the conduction in Al2O3:Cu is most likely through the interconnected Cu atoms in the network. We observed that Cu in a-Al2O3 exhibits different charge states (negative, neutral and positive) whereas the charge state of Ag in GeSe3 changes from neutral when isolated to ionic (positive) near the trapping center sites (host atoms) Chaudhuri et al. (2009).

V Acknowledgments

We thank US NSF support under grants DMR 1507670 and 1506836. We thank NVIDIA Corporation for donating a Tesla K40 GPU which was used in some of our calculations. Some of this work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1548562, using BRIDGES at the Pittsburgh Supercomputer Center under the allocation TG-DMR180083. We thank Prof. Gang Chen for valuable discussions during this work.