Structural Relaxation and Aging Scaling in the Coulomb and Bose Glass Models
We employ Monte Carlo simulations to study the relaxation properties of the two-dimensional Coulomb glass in disordered semiconductors and the three-dimensional Bose glass in type-II superconductors in the presence of extended linear defects. We investigate the effects of adding non-zero random on-site energies from different distributions on the properties of the correlation-induced Coulomb gap in the density of states (DOS) and on the non-equilibrium aging kinetics highlighted by the density autocorrelation functions. We also probe the sensitivity of the system’s equilibrium and non-equilibrium relaxation properties to instantaneous changes in the density of charge carriers in the Coulomb glass or flux lines in the Bose glass.
Disordered semiconductors near a metal-insulator transition have been the focus of experimental and theoretical investigations for the past two decades Efros1984 (); Efros1985 (); Pollak2013 (). Experimental analysis of the conductivity of a two-dimensional silicon sample revealed non-trivial non-equilibrium relaxation properties: Jaroszynski and Popović drove the system out of equilibrium by altering the gate voltage, which is equivalent to suddenly changing the charge carrier density. They monitored the conductivity relaxation and discovered that the relaxation rate of the sample depended on the waiting time popovic2007-1 (); popovic2007-2 () before performing any measurements; this served as evidence for physical aging in such systems and therefore as an indication of their complex relaxation dynamics.
Our understanding of experimental results and theoretical studies investigating such strongly correlated disordered materials stems from the Coulomb glass model that was introduced by Efros and Shklovskii Efros1984 (); Efros1985 (); Efros1975 (). This model consists of charge carriers localized at random positions and interacting via strong long-range repulsive forces. Charge carriers are exposed to pinning defect sites and they reach thermal equilibrium by rearranging themselves via variable-range hopping. The strong anticorrelations originating from the long-range repulsive interactions cause the density of states to vanish near the chemical potential, creating a soft “Coulomb gap” which in turn affects the carrier transport and equivalently the system’s conductivity. The correlation-induced soft Coulomb gap was confirmed by electron tunneling experiments in doped semiconductors Massey1995 (); Butko2000 (). In addition to the correlation-induced equilibrium properties, this model displays complex non-equilibrium relaxation dynamics that has motivated careful theoretical investigations. Monte Carlo simulations performed by Grempel et al. Grempel2004 (); Kolton2005 () discovered that the autocorrelation function measurements for different waiting times can be collapsed onto a master curve for various temperatures; an evidence of aging effects and specifically dynamical scaling Henkel2010 (). This theoretical study motivated Ovadyahu to investigate the emergence of aging effects in indium-oxide films after applying a “quench-cooling” protocol where the sample’s temperature is lowered and the system then relaxes at the new temperature for a certain waiting time Ovadyahu2006 (). The gate voltage is then abruptly changed exciting the sample and the conductance is henceforth measured to confirm that the system’s relaxation displays aging effects. Similar experimental work has been carried out by Grenet and Delahaye to study the out-of-equilibrium properties and specifically aging in disordered insulators at low temperatures, focusing on insulating granular aluminum thin films. They concluded that the sample’s relaxation indeed depends on the waiting times since the quench-cooling application Grenet2010 ().
Disordered type-II superconductors similarly require deeper understanding due to the many competing energy scales and thus rich thermodynamic phases and transport properties Blatter1994 (). The most desirable property of type-II superconductors is the zero dissipation in current flow, hence one has to prevent vortex motion due to externally applied currents to eliminate the resulting Ohmic resistance. Optimal pinning mechanisms are utilized in these materials to prevent flux flow in the presence of applied currents. Among the different pinning configurations used, columnar defects have experimentally proven to be more efficient than uncorrelated point-like disorder Civale1991 (). At low temperatures, magnetic flux lines become attached to the entire length of these linear pinning sites. This results in the emergence of the strongly pinned “Bose glass” phase Fisher-Weichman1989 (); Kwok1992 (); Nelson1992 (); Lyuksyutov1992 (); Nelson1993 () due to the localization of flux lines to these defect sites.
Physical aging has been observed in some experiments involving superconducting materials. Du et al. observed that the voltage response to an externally applied current in a superconducting 2H-NbSe sample depended on the pulse duration, and that was an evidence of the occurrence of physical aging in disordered superconducting materials Du2007 (). On the numerical front, Pleimling and Täuber employed an elastic line model and Monte Carlo simulations to study magnetic flux lines in type-II superconductors analyzing their non-equilibrium relaxation properties, starting from somewhat artificial initial conditions where straight flux lines were randomly placed in the sample Pleimling2011 (); Pleimling2015 (). The resulting complex aging features (with identical initial conditions and parameter values) were subsequently confirmed in a very different microscopic representation of the non-equilibrium vortex kinetics through Langevin molecular dynamics Dobramysl2013 (). Thereafter, Refs. Assi2015 (); Assi2016 () utilized these Langevin molecular dynamics simulations to investigate the effects of flux line density quenches on the lines’ transverse fluctuations and found that the aging scaling exponents for these transverse fluctuations depended on the choice of initial conditions Dobramysl2013 (); Assi2015 (); Assi2016 (). The structural rearrangements of these vortex lines could not be probed in the computationally accessible time scales.
The Bose glass phase in type-II superconductors and the Coulomb glass model in semiconductors share many similarities Nelson1992 (); Nelson1993 (); Dai1995 (); Tauber1995 (). At low temperatures where the magnetic vortex lines are straight and at weak magnetic fields (, with being the matching field), the Bose glass can be mapped to the Coulomb glass by considering the magnetic flux lines localized to columnar defects as particles pinned to point defects in a two-dimensional plane Nelson1993 (); Dai1995 (). The interacting flux lines in the Bose glass “hop” from their occupied defect sites to surrounding unoccupied sites. This leads to the formation of double-kinks or “superkink” excitations which are tongue-like double-kinks, each extending from one pinning site to another with a comparable energy Nelson1992 (); Nelson1993 (). Therefore, flux lines are not required to hop to nearest neighboring defect sites, which is similar to the variable-range hopping feature in the Coulomb glass. The phonon-mediated distance-depen-dent tunneling term in the Coulomb glass is replaced with an elastic energy term in the Boltzmann factor reducing the probability for long-distance double kinks in the Bose glass. Due to the remarkable similarity between the Coulomb and Bose glasses and the resulting mapping, one expects both systems to display similar properties. Analogous to the Coulomb glass, the Bose glass is composed of long-range interacting “particles” and spatial disorder, and hence this model’s density of states was confirmed to similarly display a soft gap Dai1995 (); Tauber1995 (), and the non-equilibrium relaxation dynamics display rich aging and scaling properties Shimer2010 (); Shimer2014 (). Numerical and analytical studies by Somoza et al. investigated the density of states of two-dimensional systems with logarithmic interactions at zero and finite temperatures Somoza2015 (). They confirmed that their results from the two approaches are in perfect agreement at zero temperature. Furthermore, they observed that at finite temperatures the density of states is highly dependent on the system’s temperature and can be scaled with this characteristic energy scale.
We employ Monte Carlo simulations to investigate the effects of random on-site energies and carrier density quen-ches on the density of states properties and aging scaling behavior in the Coulomb and Bose glass models. Our paper is organized as follows: The next section provides a brief theoretical background on the Coulomb glass model in disordered semiconductors and the Bose glass model in type-II superconductors in the presence of extended linear defects. It also discusses the Coulomb gap that forms in the density of states near the chemical potential due to the interaction-induced correlations. Section 3 presents an overview of the model we employ and the Monte Carlo algorithm we utilize to carry out our numerical investigations. Furthermore, it describes the protocol we perform, specifies the material parameter values we use, and presents the quantities we measure to analyze the system’s equilibrium properties and non-equilibrium and aging properties. Section 4 demonstrates the effects of adding non-zero random on-site energies on the equilibrium and non-equilibrium relaxation properties in the Coulomb and Bose glass models. Section 5 analyzes the effects of sudden changes in the density of charge carriers in the Coulomb glass model or magnetic flux lines in the Bose glass model on the density of states properties and the non-equilibrium relaxation dynamics. We finally summarize our work in the concluding section 6.
2 Theoretical Description
In this section, we describe the Coulomb glass model and the properties we analyze throughout our investigation. We thereafter move to adapt the Coulomb glass model to the Bose glass in type-II superconductors in the presence of extended linear defects.
2.1 Coulomb Glass Model
The Coulomb glass model was introduced by Efros and Shklovskii in 1974 Efros1975 () to describe localized charge carriers exposed to pinning defect centers in amorphous or doped crystalline semiconductors Efros1984 (); Efros1985 (); Efros1975 (). Since this model assumes that the localization length of charge carriers is of the same order or smaller than the mean separation between acceptor and donor sites, the carriers are confined to the randomly distributed sites. Furthermore, the carriers rearrange themselves to minimize the total interaction energy and thus reach thermal equilibrium at low temperatures. The carrier distribution rearrangements occur by variable-range hopping which is influenced by the phonon-assisted quantum tunneling between the acceptor and donor sites and thermally-induced jumps over energy barriers.
where is the electron’s charge, is the dielectric constant, and and respectively denote the position vector and occupancy of defect site ( with being the total number of defect sites), and the sites are randomly distributed on a square lattice Yu1999 (). The site occupancy is restricted to or due to the strong intra-site repulsions. represents the th site’s bare energy, and thus the first term in (1) quantifies random site energies. The second term captures the repulsive Coulomb interactions. A uniform relative charge density is introduced to maintain global charge neutrality, where is equal to the total carrier density per site or the system’s filling fraction. When we consider cases with all on-site energies , the system’s Hamiltonian (1) is characterized by particle-hole symmetry: and (where is any deviation from the half-filled system) are equivalent.
The resulting energy difference for a hop between two sites and is
where is the distance between the two sites, and the interacting site energies are given by
It is worth noting that if we replace the sites occupation numbers with the Ising spin variables , the Coulomb glass model would map into a random-site, random-field antiferromagnetic Ising model with long-range exchange interactions Davies1982 ().
2.2 Bose Glass Adaptation
The Coulomb glass model above can be adapted to describe the low-temperature properties of magnetic flux lines in type-II superconductors in the presence of strong correlated columnar disorder Nelson1993 (); Dai1995 (); Tauber1995 (). These magnetic vortices become localized at the material’s extended linear defects in the low-temperature Bose glass phase and subsequently thermal transverse fluctuations are suppressed. Therefore, the three-dimensional system becomes essentially two-dimensional by considering the cross-section of the superconductor where a pinning site is the location of the intersection of the superconductor’s “slice” transverse to the magnetic field direction and a columnar defect. The charge carriers’ role in the Coulomb glass model is played by the magnetic flux lines pinned to columnar defects in the Bose glass model.
The mutual repulsive Coulomb interaction between two occupied sites is now replaced with a modified Bessel function , which is essentially a long-range logarithmic potential screened on the scale of the London penetration length . The Bose glass effective Hamiltonian is then
where represents the system’s energy scale with the magnetic flux quantum . In the current study, we only consider a dilute low-magnetic field regime where all distances between sites are and thus .
Therefore, the energy difference for hops between two sites and in the Bose glass model is
with the interacting site energies
where (similar to the Coulomb glass model) and respectively denote the th site random on-site energy and occupation number, and is the total vortex density per site; corresponds to with being the matching field.
2.3 Coulomb Gap Properties
Efros and Shklovskii argued that the strong repulsive Cou-lomb interactions between charge carriers cause the system’s density of states to display a soft Coulomb gap in equilibrium Efros1984 (). Efros and Shklovskii’s mean-field arguments yield that the density of states vanishes at the chemical potential , which separates the low-energy filled states and the higher-energy empty states. For repulsive interactions of the form , they predict that the density of states follows a power law
for Efros1984 (); Efros1985 (); Efros1975 (); Davies1982 (); Grunewald1982 (); Xue1988 (); Grannan1993 (); Menashe2000 (); Menashe2001 (); Surer2009 (); Mobius1992 (), with a positive gap exponent . Their mean-field calculations, specifically applicable for wide random-site energy distributions Efros1979 (), predict that the gap exponent is Efros1985 (); Tauber1995 (). These mean-field arguments assume single-particle hops and that all sites are evenly distributed within the energy interval , whereas in fact multi-particle hops play an important role in the system and sites often cluster Davies1982 (). Therefore, represents a lower bound for the exponent Efros2011 (), for which Möbius, Richter, and Drittler obtained deviations from the predicted mean-field results Mobius1992 ().
The existence of the correlation-induced Coulomb gap in the density of states significantly impedes the charge carrier mobility Efros1975 () effectively lowering the sample’s conductivity: The gap exponent alters the temperature dependence of the conductivity and the system’s variable-range hopping length. In dimensions, the conductivity scales as where in the thermally-activated transport regime at low temperatures . For a constant DOS, and thus one recovers the Mott variable-range hopping exponent in three dimensions and in two dimensions. However, in the presence of the Coulomb gap, the non-zero gap exponent yields a larger : in three and two dimensions for , confirming that the Coulomb gap effectively lowers the low-temperature conductivity.
Mean-field arguments predict an exponential suppression of the density of states in the Bose glass phase, corresponding to an infinite gap exponent value. Therefore, the Coulomb gap in the Bose glass model strongly suppresses flux creep and impedes flux transport resulting in the very favorable effect of a reduced resistivity which scales as (with the driving currents ) in the thermally assisted creep regime Tauber1995 ().
3 Model and Algorithm Description
We employ Monte Carlo (MC) simulations to investigate the Coulomb glass model in disordered semiconductors and the Bose glass model in type-II superconductors with correlated columnar defects. This section specifies the MC algorithm we have utilized to study both models, the system parameters we considered, and the quantities we measured.
3.1 Model Description
Monte Carlo simulations are utilized to study the Coulomb and Bose glass models in dimensions. The simulation system is initiated by randomly placing pinning defect sites within a square lattice with periodic boundary conditions and charge carriers/flux lines that occupy sites, where the site occupancy is restricted to or . An example of the initial layout in one of our simulations is displayed in Fig. (a)a.
Then, charge carriers/flux lines attempt to hop from their corresponding occupied sites to unoccupied sites. The transition rate from site to site is determined by
where denotes a microscopic time scale, represents the spatial extension of the charge carrier wave function or thermal wandering of the magnetic flux lines, and given by Eq. (2) for the Coulomb glass model or (5) for the Bose glass model is the difference between the energies of the two sites in question. Therefore, the success of hops is determined by two factors in (8): an exponential term capturing the strongly distance-dependent phonon-mediated tunneling in semiconductors and vortex superkink proliferation in type-II superconductors, and a Metropolis term describing thermally-activated jumps through energy barriers. These hops allow the system in Fig. (a)a to relax to a new configuration shown in Fig. (b)b after MC time steps (the MC algorithm is outlined in the following section).
3.2 Monte Carlo Algorithm
One Monte Carlo time step (MCS) in our simulations occurs as follows:
Choose a random occupied defect site .
Choose a random unoccupied defect site with probability .
Pick a random number , and attempt a charge carrier hopping from to . The probability of the transition success is given by the Metropolis term, , i.e., move the charge carrier from defect site to site if . If the hopping attempt fails, return to the first step.
One MCS consists of N iterations of the outlined steps.
3.3 Simulation Protocol
The first system we study as a continuation to the work in Refs. Shimer2010 (); Shimer2014 () is the initial high-temperature quench with completely random initial conditions. The simulation lattice with the corresponding number of defect sites and carriers is prepared, and the system relaxes for a sufficiently long initial relaxation time for both models.
The second system we consider involves an abrupt change in the number of charge carriers, where carriers are randomly added or removed from the sample. This serves as an attempt to bridge our computational work with experimental research and to better understand the system’s sensitivity to the choice of initial conditions. For this part, we let the system relax for a long time beyond any microscopic time scale . At MCS, the number of charge carriers is instantaneously increased from to or lowered from to . Then, we let the resulting system relax up to a waiting time measured after the quench. A snapshot of the system is taken at and the two-time carrier density correlation function defined in Section 3.5 is measured at later times . A similar study of flux line density quenches is performed to analyze structural rearrangements of flux lines in the Bose glass phase in disordered type-II superconductors.
3.4 System Parameters
All simulation distances are measured in units of the average distance between neighboring defect sites . Simulation energies in the Coulomb glass model are measured in units of the typical Coulomb energy , whereas energies in the Bose glass model are in units of . Furthermore, all simulation times are measured in units of Monte Carlo time steps (MCS). The spatial extension of the charge carrier wave function is chosen to be Grempel2004 (); Kolton2005 () following the careful investigation of different values of in Ref. Shimer2010 (). In the Bose glass model, the London penetration length is set to .
Earlier work by Efros et al. Efros2011 () focused on random site energies from a uniform distribution of certain width, while here we investigate different energy distributions to conclude their effects on the system’s properties. We consider cases with zero random on-site energies () and others where the on-site energies are taken from either a normalized Gaussian or a flat distribution of zero mean and different widths .
All other simulation parameters are chosen in accordance with the findings in Ref. Shimer2010 (). The simulation lattice length is , since the study of various lengths proved the absence of measurable finite-size effects. The system’s ambient temperature is set to , since the temperature range that yields reasonable dynamics is : the system’s kinetics slowed down too much below , and no aging was observed above . The filling fraction is kept within the range since the dynamics in systems with filling fractions or freezes out in the accessible simulation time scales. Periodic boundary conditions are utilized, and each simulation configuration shown was analyzed with to independent simulation runs.
3.5 Measured Quantities
The density of states (DOS) , with the interacting site energies given by (3) for the Coulomb glass model and (6) for the Bose glass model, is measured by finding the number of sites that are located in each energy bin, with the bin size chosen as . We also study the soft gap that forms in the DOS and the speed of its formation over time given different simulation configurations. Mean-field calculations predict the DOS to display a power law near the chemical potential , which directs us to measure the gap exponent and compare it to the mean-field predicted .
The non-equilibrium relaxation dynamics and response to the high-temperature quench and the charge carrier (or magnetic flux line) density abrupt changes are analyzed using the normalized two-time carrier density autocorrelation function
where is the occupation number of defect site at time , is the waiting time, and is the filling fraction. Furthermore, denotes an average over all sites, followed by an average over several thousand independent disorder realizations and thermal noise.
4 Random On-site Energy Effects
Via Monte Carlo simulations, we have studied the emergence of the soft gap in the density of states in the Coulomb and Bose glass models in two dimensions. We compare the gap’s properties and non-equilibrium relaxation dynamics of systems without random on-site energies to those with on-site energies from different distributions. Section 4.1 discusses the effects of adding non-zero on-site energies on the density of states and the Coulomb gap’s characteristics in the Coulomb and Bose glass models, while Section 4.2 highlights the ensuing non-equilibrium relaxation dynamics and resulting aging scaling exponents. In this section, there are random pinning sites available to charge carriers/flux lines.
4.1 Coulomb Gap Properties
4.1.1 Coulomb Glass Model
We first investigate the Coulomb gap formation in the Coulomb glass model under random initial conditions in the absence of random on-site energies, where a fixed number of charge carriers is exposed to randomly distributed defect sites. At the beginning of our simulations before any hops occur, the density of states has a maximum at the chemical potential . At MCS, a local minimum in the DOS already starts to form as a clear indication for the occupation of some low-energy states due to the charge carriers hopping between sites. The DOS in Fig. 2 is almost totally suppressed at the chemical potential and the Coulomb gap is pronounced. This confirms the occupation of the low-energy states by the present charge carriers, which is pictured by the peak located at , and the separation from the unoccupied higher-energy states at .
We have monitored the evolution of the DOS in time in the absence of on-site energies to gauge the effect of the inclusion of such energies on the speed of the gap’s formation. To this end, we studied systems without on-site energies and others with random on-site energies chosen from either a normalized Gaussian or a flat distribution of different widths. Fig. 3 shows the gap in the DOS starting to form within MCS, and the value of decreasing as time progresses until (an almost) total suppression is observed in our two-dimensional Coulomb glass model. The three-dimensional system was observed in Refs. Shimer2010 (); Shimer2014 () to display a faster total suppression of the Coulomb gap than the two-dimensional system.
Fig. 3 confirms that the inclusion of random on-site energies from a Gaussian distribution of different widths (here and are utilized) does not influence the speed of the gap’s formation and it does not affect the suppression of the DOS at the chemical potential at long times. On the other hand, one observes that on-site energies from a flat distribution of zero mean and width might be slowing down the suppression of the DOS at , where reaches a plateau value after MCS of our accessible time scales as seen in Fig. 3. It is possible that the inclusion of this energy distribution might have caused the density of states to freeze out in our small-sized Coulomb glass system for the accessible time scales considered here.
Comparing Fig. 4 to Fig. 2 without on-site energies at the same time ( MCS) asserts that on-site energies influence the arrangement of charge carriers and their pinning behavior to defect sites, thereby eliminating the symmetry of the DOS around for half-filling and more importantly delaying the time when total suppression of the DOS at is reached.
To better understand the effect of non-zero random on-site energies on the Coulomb gap, we have studied the asymptotic behavior of the density of states near the chemical potential. As discussed earlier, the DOS near the chemical potential diminishes via a power law and mean-field arguments compute the gap exponent as for the Coulomb repulsive interactions of the form . This implies that the mean-field gap exponent for the Coulomb glass model in two dimensions is . We study the asymptotic behavior of the density of states near the chemical potential for the cases characterized by the absence and later presence of different on-site energies and observe that a power law behavior is displayed in the DOS near ; an example is shown in Fig. 5.
As indicated in Fig. 5, the computed effective gap exponent for the system without random on-site energies is found to be , which markedly differs from the mean-field prediction . As discussed in Section 2, it is expected for the exponent computed in these simulations to be different from that predicted by mean-field arguments. Furthermore, the mean-field prediction pertains to systems with strong disorder (wide random site energy distributions) Efros1979 (). We have measured the gap exponent in cases with different on-site energy distributions, and the values we found are summarized in Table 1.
Note that in order to estimate the characteristic error associated with the slope of the linear fit of our data, thus accounting for varying fluctuation sizes in different portions of the data, we divide the dataset into multiple regions. We first use a simple linear regression to obtain a line of best fit for the entire data range, producing the major line. We then divide the dataset into equal-sized regions along the x-axis and regressively fit a straight line to each region and note the slope of the line and standard deviation in the slope, two numbers that define an interval for the slope. This gives us two slope values associated with each region, i.e., the upper and lower bound of the interval. We then compare all the slope values obtained to the slope of the major line and report the largest difference as the characteristic error associated with the linear fit to our data.
The flat random on-site energy distribution has a much more pronounced influence on the value of the computed exponent than the Gaussian energy distributions. Including random on-site energies from this flat distribution decreases the gap exponent and renders it closer to the mean-field predictions, as expected Efros1979 (). Sharp energy distributions are less effective due to the interaction-induced correlations, and one obtains results that are very similar to the case when (we confirmed this with flat and Gaussian distributions of zero mean and width ). In contrast, broader distributions with variances on the same scale as the correlation-induced Coulomb gap’s width have a more pronounced effect on the density of states around the chemical potential. A lower effective gap exponent than the case without on-site energies implies a lower , which in turn translates to an increased conductivity and thus charge carrier transport since conductivity scales as .
|(a)||(b) Gaussian||(c) Gaussian||(d) Flat|
4.1.2 Bose Glass Model
We have similarly monitored the density of states in the Bose glass model in two dimensions and analyzed the Coulomb gap’s properties that characterize this system. A very pronounced difference with the Coulomb glass model is that the soft gap forms quicker and the DOS is suppressed faster, i.e., in shorter time, in about MCS. The total suppression of the DOS at for the Bose glass can be observed in Fig. 6.
The Coulomb gap forms much faster and is broader than that displayed by systems with Coulomb repulsive interactions. Analyzing the evolution of in time confirms this observation and the direct effect of the addition of random on-site energies from various distributions. Comparing the speed of formation of the Coulomb gap in systems with logarithmic interaction displayed in Fig. 7 to that in systems with Coulomb interactions in Fig. 3 proves that the logarithmic repulsive interactions drive the gap to form much faster.
Aside from the faster formation of the gap, the DOS is also strongly suppressed for systems without on-site energies and energies from flat and Gaussian distributions alike. Therefore, one concludes that the inclusion of on-site energies from different distributions does not change the overall properties and speed of formation of the gap in the Bose glass model, unlike our findings in the two-dimensional Coulomb glass model where energies from a flat distribution affect the suppression of . To further investigate the effects of random on-site energies on the density of states in the Bose glass model, we have checked for the power law fit and the gap exponent in this case. Similar to the Coulomb glass model in Section 4.1.1, we observe that follows a power law in the vicinity of the chemical potential .
However, the gap exponents computed in the Bose glass are large when compared to those in the Coulomb glass. One has to keep in mind that the suppressed motion within the Bose glass model leads to a slow relaxation and thus a large gap exponent. Calculating the gap exponent from the power law best fit allows us to study the influence of random on-site energies on the density of states. We have computed the Coulomb gap exponents from the power law best fit in the Bose glass system in the absence of any on-site energies and in the presence of random energies from a flat distribution of width and a Gaussian distribution of widths or . The summary of our findings is presented in Table 2.
|(a)||(b) Gaussian||(c) Gaussian||(d) Flat|
On-site energies from a Gaussian distribution of different widths show a much weaker effect on the gap exponent than the flat distribution, similar to our findings in the Coulomb glass model. The introduction of random on-site energies from the utilized flat distribution increases the effective gap exponent rendering it closer to the infinite value that mean-field arguments predict in the Bose glass phase.
4.2 Non-Equilibrium Relaxation Dynamics
In addition to our investigations of the effects of random on-site energies on the density of states, the Coulomb gap formation, and the effective gap exponent, we have analyzed the corresponding effects on the non-equilibrium relaxation properties of the systems studied following a quench from a completely uncorrelated, high-temperature initial state. Section 4.2.1 presents our findings for the Coulomb glass model in disordered semiconductors, and Section 4.2.2 summarizes our results for the Bose glass model in type-II superconductors in the presence of extended linear defects.
4.2.1 Coulomb Glass Model
We have investigated the non-equilibrium relaxation dynamics in the Coulomb glass model by measuring the two-time carrier density autocorrelation function and its properties.
In the absence of random on-site energies, the system shows slow dynamics and breaking of time translation invariance due to its dependence on the waiting time, as seen in Fig. 8. The slow relaxation of this system with Coulomb repulsive interactions is apparent in Fig. 8 since the system has not reached a steady state even after MCS. The system displays two different relaxation regimes, similar to those shown in Fig. 9 for the system with Gaussian on-site energies.
The first regime is characterized by a fast relaxation produced by the displacements of charge carriers around their initial positions and the carrier clusters’ group motion. On the other hand, the second regime displays a slower relaxation that is dominated by the divisions of charge carrier clusters and the resulting relaxation of the system towards a steady state. The slow dynamics and breaking of time translation invariance are two characteristics of physical aging, while the third property to be considered here is dynamical scaling Henkel2010 ().
The carrier density autocorrelation function in the aging scaling regime obeys the general scaling form
where is an aging scaling exponent. When in the case of full aging, the scaling function is a power law , where is the autocorrelation exponent and is the dynamic exponent Henkel2010 ().
Investigating the dynamical scaling property in this regime, we observe that scaling collapse is obtained with the exponents and for the two-dimensional Coulomb glass model without on-site energies, see Fig. (a)a. The goal of this section is to compute these scaling exponents after the addition of on-site energies from different distributions and compare them to the above case. We focus on energies from flat and Gaussian distributions, and we present the measured scaling exponents in Table 3.
|(a)||(b) Gaussian||(c) Gaussian||(d) Flat|
From these results, we conclude that the addition of random on-site energies from a Gaussian distribution of different widths does not noticeably affect the aging scaling exponents. On the other hand, on-site energies from a flat distribution of width highly influence both exponents, see Fig. (b)b. The introduction of broad on-site energy distributions of widths on the same scale as the correlation-induced Coulomb gap’s width drives the system to display a faster relaxation dynamics, as confirmed by the higher .
4.2.2 Bose Glass Model
A similar study of the non-equilibrium relaxation properties of the Bose glass model in two dimensions has been performed, where we have analyzed the effects of random on-site energies on the density autocorrelation function and its scaling form and exponents.
Starting with the Bose glass model with zero on-site energies, slow dynamics is also a signature of such a system and the autocorrelation function displays a two-regime relaxation similar to that of the Coulomb glass model in Fig. 9. Investigating the universality of the scaling exponents, we study the full aging scaling forms for the various on-site energy distributions.
To monitor the effects of on-site energies on the two-dimensional Bose glass system’s non-equilibrium relaxation properties, we similarly investigate any discrepancies present between the scaling exponents in the case with zero on-site energies and those in the cases with different energy distributions. The results are summarized in Table 4.
|(a)||(b) Gaussian||(c) Gaussian||(d) Flat|
Our results in Table 4 confirm that adding random on-site energies from a Gaussian distribution of different widths has a minute effect on the scaling exponents and , see Fig 11, which implies that this addition hardly modifies the relaxation dynamics of the Bose glass system with essentially logarithmic interactions. On the other hand, it is evident that a flat on-site energy distribution of width has a pronounced effect on the non-equilibrium relaxation properties of this system. Both aging scaling exponents and assume considerably larger values as a response to the introduction of non-zero random on-site energies from a flat distribution, which hence imposes a faster relaxation on the system analogous to our results for the Coulomb glass model.
5 Density Quench Effects
Previous studies addressed the Coulomb gap and non-equilibrium relaxation properties under random initial conditions Grempel2004 (); Kolton2005 (); Shimer2010 (); Shimer2014 (). It is worthwhile to analyze the effect of abrupt changes in the density of charge carriers in the Coulomb glass model or analogously flux lines in type-II superconductors with linear defects on the Coulomb gap’s properties and aging scaling exponents to probe the system’s sensitivity to sudden changes in parameters. These quenches are experimentally achieved by quickly switching the gate voltage for semiconductors or the external magnetic field in type-II superconductors. We have utilized Monte Carlo simulations to study the effects of density quenches on the Coulomb gap in Section 5.1 and on the aging scaling behavior and exponents in Section 5.2, all in the absence of random on-site energies ().
5.1 Coulomb Gap Properties
5.1.1 Coulomb Glass Model
In the Coulomb glass model in two dimensions, we investigate the effects of suddenly increasing the filling fraction from to (up-quench) or alternatively decreasing the filling fraction from to (down-quench) on the Coulomb gap and its evolution in time.
We first keep the density of charge carriers fixed with the filling fraction and let the system relax in time for MCS until the Coulomb gap forms at the chemical potential where for half-filled systems (this is when the simulation clock is reset to ). After confirming that the density of states displays a relaxed Coulomb gap, we perform the density quench. Discussing first the case of the density up-quench and measuring the DOS, we have to keep in mind that the system is not at half-filling anymore but at , which implies that we can expect a shift in the DOS as compared to the case with a fixed density. To test this expectation, we follow the DOS evolution in time.
At the moment of the charge carrier density up-quench, Fig. (a)a shows that the peak on the left denoting the lower-energy states reduces while the opposite peak increases. This occurs since new charge carriers are abruptly introduced to the system which in turn increases the typical Coulomb repulsion, and thus an enhanced peak on the higher-energy states initially emerges. As the system relaxes away from this initial arrangement and towards equilibrium, the new charge carriers start exploring their surroundings for energetically more favorable pinning sites. This results in a change of the shape of the DOS through a flip in the asymmetry, compare Fig. (a)a with Fig. (b)b.
Since after the density quench, this disappearance of the symmetry that characterizes the half-filled system is expected. In the case of an abrupt increase in away from half-filling, the equilibrium chemical potential is shifted to a new value .
The case of density down-quench is similarly performed, where the filling fraction is abruptly lowered from to . At the moment of the quench at , the peak at the lower energies is enhanced, which is expected due to the reduction in the number of charge carriers which implies a decrease in the typical Coulomb repulsion. A similar but reversed asymmetry to that in Fig. (b)b therefore appears at MCS after the down-quench. As the system relaxes, a smaller number of charge carriers occupies the original pinning defect sites, and hence the asymmetry is reversed to show an increased peak for the unoccupied states.
Another property that can be utilized to investigate the effects of instantaneous density changes is the power law fit for the density of states in the vicinity of the chemical potential. We compute the effective gap exponent at the same time for the cases when the filling fraction is kept fixed at half or suddenly raised/lowered away from half-filling.
|(a) Fixed density||(b) Up-quench||(c) Down-quench|
From Table 5, one obtains the effective gap exponents for the three systems with different final filling fractions. It is worth noting that since the Coulomb gap exponent is an equilibrium property, one expects to obtain similar results when either starting from random initial conditions with filling fraction or abruptly changing the filling fraction from to and thereafter letting the system relax towards equilibrium.
5.1.2 Bose Glass Model
Having established that the DOS configuration is affected by sudden changes in the charge carrier density in the Coulomb glass model, we aim to carry out the same investigation for the Bose glass model with logarithmic interactions. This study is significant because it analyzes the effects of abrupt changes in the external magnetic field on the spatial rearrangements of magnetic vortex lines in type-II superconductors, while our earlier work in Ref. Assi2015 (); Assi2016 () focused on the flux lines height fluctuations and the structural relaxation time regime remained inaccessible.
Similar to the study above, we start with the Bose glass model with a fixed filling fraction , and we obtain the density of states in Fig. 6, where a Coulomb gap appears and the DOS is broad and totally suppressed at the chemical potential , which implies that the system here is relaxed beyond microscopic time scales. In the case of suddenly increasing the system’s filling fraction from to , one notices the same behavior in the density of states in the Bose glass where the system starts with an enhanced peak at the higher-energy states, similar to Fig. (a)a. Hence, the symmetry that is witnessed with a fixed number of magnetic flux lines at half-filling is broken due to the initial introduction of new magnetic flux lines into the system with a constant number of pinning centers. As the system relaxes to MCS, the newly-added flux lines have explored the surrounding landscape and become pinned to the defect sites, and hence a peak on occupied sites is now enhanced. When we consider the consequences of abruptly reducing the number of flux lines due to changing the filling fraction from to , we observe the same reversed asymmetry that is witnessed in the down-quench case in the Coulomb glass system.
|(a) Fixed density||(b) Up-quench||(c) Down-quench|
We studied the effects sudden changes in the number of flux lines have on the effective gap exponent in the Bose glass model. The results are summarized in Table 6. The effective gap exponents for the three different final filling fractions in Table 6 are almost equal, which implies that this equilibrium property is similar in systems with the studied filling fractions , , and in the Bose glass model.
5.2 Non-Equilibrium Relaxation Dynamics
In addition to studying the effects of charge carrier/flux line density quenches on the density of states and the soft Coulomb gap, we analyzed the effects of sudden addition/removal of carriers on the non-equilibrium relaxation properties and the aging scaling exponents: Section 5.2.1 presents the results of this study in the Coulomb glass model and Section 5.2.2 summarizes our findings for the Bose glass model, all in the absence of random on-site energies ().
5.2.1 Coulomb Glass Model
In the two-dimensional Coulomb glass model, we relax the system for MCS (this is when the simulation clock is reset to ) to then quench the density of charge carriers and measure the resulting two-time carrier density autocorrelation function for different waiting times measured after the quench.
In the case where the density is kept fixed at , the dynamics in Fig. 13 becomes less influenced by the choice of waiting times due to the long initial relaxation time that this finite system undergoes. The density autocorrelation function of this relaxed system with a fixed does not acquire the dynamical scaling property, regardless of the choice of scaling exponents, possibly due to the finite-size effects that emerge at the time scales considered. When the charge carrier density suddenly changes, we observe different features to the case with a fixed filling fraction .
For both cases when (a) the charge carrier density abruptly increases from to and (b) the density abruptly decreases from to , the dynamics is again highly dependent on the waiting times chosen as can be seen in Fig. (a)a. Furthermore, the case of density down-quench shows a similar relaxation dynamics to that of the up-quench case in Fig. (a)a.
Dynamical scaling is not present in the relaxed system with a fixed carrier density, while a similar system with random initial conditions displays dynamical scaling with the scaling exponents and . Investigating the dynamical scaling property after the influence of sudden changes in the density of charge carriers reveals that aging scaling is again a property of the system. Due to the similar relaxational dynamics in the cases of density up-quench and down-quench, we confirm that the dynamical aging scaling exponents and are equal in both cases. This dynamical scaling behavior is displayed in Fig. (b)b, from which we obtain the values and .
5.2.2 Bose Glass Model
The relaxational dynamics of the Bose glass model in two dimensions was also analyzed with a similar goal of studying the effects of sudden changes in the density of magnetic flux lines.
When the density of flux lines is kept fixed at and the system is relaxed for a sufficiently long initial time MCS, the relaxation of the autocorrelation function is very weakly dependent on the waiting times chosen, see Fig. 15.
On the other hand, when the flux line density is suddenly increased from to or decreased from to , we confirm that the system is again dependent on the waiting times . Moreover, both systems similarly display an initial fast relaxation regime and a slower later regime, see Fig. 16.
Dynamical scaling is not present in the relaxed system with a fixed flux line density, while a similar system with random initial conditions displays dynamical scaling with the scaling exponents and . In situations with abrupt increases or decreases in the magnetic flux line density, the following values for these scaling exponents are found: and .
The equal aging scaling exponents and in the Bose glass systems when the magnetic flux line density is abruptly increased or decreased implies that this system undergoes similar relaxational dynamics when new flux lines are added to the sample or after some flux lines are removed from the sample due to particle-hole symmetry. This observation was also made in the Coulomb glass model with sudden changes in the charge carrier density.
We have employed the Coulomb glass model in disordered semiconductors and adapted it to the Bose glass model in type-II superconductors in the presence of extended linear defects to investigate the density of states and the non-equilibrium relaxation properties in the two-dimensional Coulomb and Bose glass systems via Monte Carlo simulations.
Earlier investigations focused on these systems with zero on-site energies Grempel2004 (); Kolton2005 (); Shimer2010 (); Shimer2014 (), thus one goal of this study is to analyze the effects of adding random on-site energies from different distributions into the Coulomb and Bose glass systems with a fixed charge carrier/flux line density. We conclude that adding on-site energies from a Gaussian distribution of zero mean and of various widths does not noticeably affect the speed of formation of the Coulomb gap or the effective gap exponent in both models. On the other hand, a flat distribution of a similar width to the one utilized for the Gaussian distribution causes the suppression of the density of states at the chemical potential to become slower in the Coulomb glass model (keeping in mind the possibility that the flat on-site energy distribution might have caused the density of states to freeze out in our small-sized Coulomb glass system in the accessible simulation time scales), while not affecting this feature in the Bose glass model. Furthermore, adding random on-site energies from this flat distribution causes the gap exponent’s effective value to approach the mean-field predictions in both models. Such broad on-site energy distributions of widths on the same scale as the correlation-induced Coulomb gap width affect the density of states near the chemical potential in these two models.
The non-equilibrium relaxation dynamics shows the effects of non-zero random on-site energies as well: The Gaussian-distributed on-site energies do not change the aging scaling exponents, while a wide flat distribution has a pronounced effect on these exponents in both systems. Including random on-site energies from a flat distribution into the Coulomb and Bose glass models causes their relaxation to become faster. This change in the aging scaling exponents implies that these exponents are not universal and in fact depend on various microscopic system parameters Shimer2010 (); Shimer2014 ().
Another goal of this study is to analyze the effects of abrupt changes in the density of charge carriers in the Coulomb glass model or magnetic flux lines in the Bose glass model to gauge the systems’ response to sudden perturbations. In the cases of density up-quench and down-quench in both models, the system first responds to the quench by shifting the whole DOS curve to higher (up-quench) or lower (down-quench) energies and displaying an asymmetry in the density of states around its minimum value. When the modified system reaches equilibrium, the larger peak shifts to the opposite side of the new chemical potential due to the relaxation of the collective system with the newly-added carriers or the relaxation of the same system after removing some carriers. Since the Coulomb gap exponent is an equilibrium property, one should of course obtain the same exponent in a system with a specific filling fraction regardless of the utilized initial conditions.
We analyzed the effects of density quenches on the non-equilibrium relaxation dynamics and the aging properties in the Coulomb and Bose glass models. Systems with a fixed density relax for a long time until time translation invariance is displayed by the carrier density autocorrelation function. We observe that dynamical scaling is not a property of finite Coulomb and Bose glass systems at very long times when the density is fixed. On the other hand, when the density is suddenly raised or lowered, the system displays dynamical scaling and aging scaling exponents that are equal in the cases of density up-quench and down-quench in both models. The equal aging scaling exponents when the density of charge carriers or flux lines abruptly changes confirms that the corresponding systems undergo similar relaxational dynamics as a response to the addition of new carriers or the removal of an equivalent number of carriers due to particle-hole symmetry. Within error bars, the aging scaling exponents in the Coulomb and Bose glass models with random initial conditions and after sudden changes in the carrier density coincide, suggesting that structural rearrangements of charge carriers/flux lines are governed by the low-energy equilibrium site energy distributions inside the Coulomb gap.
We gladly acknowledge stimulating and helpful discussions with Boris Shklovskii, Matt Shimer, and Shadi Esmaeili. This work was supported by the U.S. Department of Energy, Office of Basic Energy Sciences, Division of Materials Sciences and Engineering, under Grant No. DE-FG02-09ER46613.
Authors Contribution Statement
Hiba Assi carried out the numerical work, analyzed the data, and wrote the first manuscript draft. Harshwardhan Chaturvedi assisted in the computational part and coding. Michel Pleimling co-advised the research, contributed to data analysis and physical discussions, and assisted in finalizing the paper. Uwe C. Täuber proposed the research topic, supervised the research progress and data analysis, and edited the first manuscript draft.
- (1) B. I. Shklovskii and A. L. Efros, Electronic Properties of Doped Semiconductors (Springer, New York, 1984).
- (2) A. L. Efros and M. Pollak, Electron-Electron Interactions in Disordered Systems (North Holland, Amsterdam, 1985).
- (3) M. Pollak, M. Ortuño, and A. Frydman, The Electron Glass (Cambridge University Press, New York, 2013).
- (4) J. Jaroszynski and D. Popović, Phys. Rev. Lett. 99, 046405 (2007).
- (5) J. Jaroszynski and D. Popović, Phys. Rev. Lett. 99, 216401 (2007).
- (6) A. L. Efros and B. I. Shklovskii, J. Phys. C: Solid State Phys. 8, L49 (1975).
- (7) J. G. Massey and M. Lee, Phys. Rev. Lett. 75, 4266 (1995).
- (8) V. Yu. Butko, J. F. DiTusa, and P. W. Adams, Phys. Rev. Lett. 84, 1543 (2000).
- (9) D. R. Grempel, Europhys. Lett. 66, 854 (2004).
- (10) A. B. Kolton, D. R. Grempel, and D. Domínguez, Phys. Rev. B 71, 024206 (2005).
- (11) M. Henkel and M. Pleimling, Nonequilibrium Phase Transitions, Volume 2 – Ageing and Dynamical Scaling far from Equilibrium (Springer, Heidelberg, 2010).
- (12) Z. Ovadyahu, Phys. Rev. B 73, 214204 (2006).
- (13) T. Grenet and J. Delahaye, Eur. Phys. J. B 76, 229 (2010).
- (14) For a now classic review, see: G. Blatter, M. V. Feigel’man, V. B. Geshkenbein, A. I. Larkin, and V. M. Vinokur, Rev. Mod. Phys. 66, 1125 (1994).
- (15) L. Civale, A. D. Marwick, T. K. Worthington, M. A. Kirk, J. R. Thompson, L. Krusin-Elbaum, Y. Sun, J. R. Clem, and F. Holtzberg, Phys. Rev. Lett. 67, 648 (1991).
- (16) M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. Fisher, Phys. Rev. B 40, 546 (1989).
- (17) W. K. Kwok, S. Fleshler, U. Welp, V. M. Vinokur, J. Downey, G. W. Crabtree, and M. M. Miller, Phys. Rev. Lett 69, 3370 (1992).
- (18) D. R. Nelson and V. M. Vinokur, Phys. Rev. Lett. 68, 2398 (1992).
- (19) I. F. Lyuksyutov, Europhys. Lett. 200, 273 (1992).
- (20) D. R. Nelson and V. M. Vinokur, Phys. Rev. B 48, 13060 (1993).
- (21) X. Du, G. Li, E. Y. Andrei, M. Greenblatt, and P. Shuk, Nature Phys. 3, 111 (2007).
- (22) M. Pleimling and U. C. Täuber, Phys. Rev. B 84, 174509 (2011).
- (23) M. Pleimling and U. C. Täuber, J. Stat. Mech. (2015) P09010.
- (24) U. Dobramysl, H. Assi, M. Pleimling, and U. C. Täuber, Eur. Phys. J. B 86, 228 (2013).
- (25) H. Assi, H. Chaturvedi, U. Dobramysl, M. Pleimling, and U. C. Täuber, Phys. Rev. E 92, 052124 (2015).
- (26) H. Assi, H. Chaturvedi, U. Dobramysl, M. Pleimling, U. C. Täuber, to appear in Molecular Simulation (2016) [arXiv:1509.02227].
- (27) U. C. Täuber, H. Dai, D. R. Nelson, and C. M. Lieber, Phys. Rev. Lett. 74, 5132 (1995).
- (28) U. C. Täuber and D. R. Nelson, Phys. Rev. B 52, 16106 (1995).
- (29) M. T. Shimer, U. C. Täuber, and M. Pleimling, Europhys. Lett. 91, 67005 (2010).
- (30) M. T. Shimer, U. C. Täuber, and M. Pleimling, Phys. Rev. E 90, 032111 (2014).
- (31) A. M. Somoza, M. Ortuño, T. I. Baturina, and V. M. Vinokur, Phys. Rev. B 92, 064201 (2015).
- (32) C. C. Yu, Phys. Rev. Lett. 82, 4074 (1999).
- (33) J. H. Davies, P. A. Lee, and T. M. Rice, Phys. Rev. Lett. 49, 758 (1982).
- (34) M. Grünewald, B. Pohlmann, L. Schweitzer, and D. Würtz, J. Phys. C: Solid State Phys. 15, L1153 (1982).
- (35) W. Xue and P. A. Lee, Phys. Rev. B 38, 9093 (1988).
- (36) E. R. Grannan and C. C. Yu, Phys. Rev. Lett. 71, 3335 (1993).
- (37) D. Menashe, O. Biham, B. D. Laikhtman, and A. L. Efros, Europhys. Lett. 52, 94 (2000).
- (38) D. Menashe, O. Biham, B. D. Laikhtman, and A. L. Efros, Phys. Rev. B 64, 115209 (2001).
- (39) B. Surer, H. G. Katzgraber, G. T. Zimanyi, B. A. Allgood, and G. Blatter, Phys. Rev. Lett. 102, 067205 (2009).
- (40) S. D. Baranovskii, A. L. Efros, B. L. Gelmont, and B. I. Shklovskii, J. Phys. C: Solid State Phys. 12, 1023 (1979).
- (41) A. Möbius, M. Richter, and B. Drittler, Phys. Rev. B 45, 11568 (1992).
- (42) A. L. Efros, B. Skinner, and B. I. Shklovskii, Phys. Rev. B 84, 064204 (2011).