text:
# Chemical Potential Shift in Nd2-xCexCuO4: Contrasting Behaviors of the Electron- and Hole-Doped Cuprates

## Abstract

We have studied the chemical potential shift in the electron-doped superconductor Nd2-xCexCuO4 by precise measurements of core-level photoemission spectra. The result shows that the chemical potential monotonously increases with electron doping, quite differently from La2-xSrxCuO4, where the shift is suppressed in the underdoped region. If the suppression of the shift in La2-xSrxCuO4 is attributed to strong stripe fluctuations, the monotonous increase of the chemical potential is consistent with the absence of stripe fluctuations in Nd2-xCexCuO4. The chemical potential jump between Nd2CuO4 and La2CuO4 is found to be much smaller than the optical band gaps.

High-$`T_c`$ superconductivity occurs when the parent antiferromagnetic (AF) insulator with the CuO2 plane is doped with holes or electrons. In the p-type materials, the long-range AF order vanishes for a slight amount of hole doping whereas in the n-type materials, the AF order persists up to a high doping concentration of $``$0.14 electrons per Cu and the superconducting (SC) doping range is much narrower . The p-type materials show $`T`$-linear in-plane electrical resistivity  and split neutron peaks around $`𝐪=(\pi ,\pi )`$ indicating incommensurate spin fluctuations  whereas the n-type materials show $`T^2`$ dependence of the in-plane resistivity  and $`(\pi ,\pi )`$ commensurate spin fluctuations . In order to elucidate the mechanism of high-$`T_c`$ superconductivity, it is very important to clarify the origin of the similarities and the differences between the p-type and the n-type materials.

In this Letter, we report on a study of the chemical potential shift in Nd2-xCexCuO4 (NCCO) as a function of doped electron concentration. The shift can be deduced from the core-level shifts in photoemission spectra because the binding energy of each core level is measured relative to the chemical potential $`\mu `$. In a previous study , we found that in La2-xSrxCuO4 (LSCO) the chemical potential shift is unusually suppressed in the underdoped region and attributed this observation to the strong stripe fluctuations which exist in this system. As for the chemical potential jump between La2CuO4 and Nd2CuO4, which would represent the band gap of the parent insulator, it was estimated to be at most 300 meV in previous valence-band photoemission studies , which is much smaller than the 1.5–2.0 eV charge-transfer (CT) gap of the parent insulator estimated from optical studies .

High-quality single crystals of NCCO ($`x=`$ 0, 0.05, 0.125 and 0.15) were grown by the traveling-solvent floating-zone method as described elsewhere . Uncertainties in the Ce concentration were $`\pm 0.01`$. For $`x=0.15`$, both as-grown and reduced samples were measured while for the other compositions only as-grown samples were measured. The as-grown samples were all antiferromagnetic and did not show superconductivity. Only the $`x=0.15`$ sample showed superconductivity after reduction in an Ar atmosphere and its $`T_c`$ was $`25`$ K.

X-ray photoemission spectroscopy (XPS) measurements were performed using both the Mg $`K\alpha `$ ($`h\nu =1253.6`$eV) and Al $`K\alpha `$ ($`h\nu =1486.6`$eV) lines and a hemispherical analyzer. All the spectra were taken at liquid-nitrogen temperature ($`80`$ K) within $`40`$ minutes after scraping. We did not observe a shoulder on the higher binding energy side of each O 1s peak, indicating the high quality of the sample surfaces free from degradation. Although the energy resolution was about $`0.8`$ eV for both $`K\alpha `$ lines, we could determine the core-level shifts with an accuracy of about $`\pm 50`$ meV because most of the spectral line shapes did not change with $`x`$. In XPS measurements, a high voltage of $`>1`$ kV is used to decelerate photoelectrons, and it is usually difficult to stabilize the high voltage with the accuracy of $`100`$ meV. In order to overcome this difficulty, we directly monitored the voltage applied to the outer hemisphere and the retarding fringe, and confirmed that the uncertainty could be reduced to less than $`10`$ meV. To eliminate other unexpected causes of errors, we measured the $`x=0.05`$ sample as a reference just after the measurement of each sample.

Figure 1 shows the XPS spectra of the O 1s, Nd 3$`d_{5/2}`$ and Cu 2$`p_{3/2}`$ core levels taken with the Al $`K\alpha `$ line. Here, the integral background has been subtracted and the intensity has been normalized to the peak height . The Nd 3$`d_{5/2}`$ spectra are composed of the 3$`d_{5/2}`$4$`f^4`$ and 3$`d_{5/2}`$4$`f^3`$ final-state components, where  denotes a ligand hole, and O KLL Auger signals overlap them. The Cu 2$`p_{3/2}`$ spectra are composed of the 2$`p_{3/2}`$3$`d^{10}`$ and 2$`p_{3/2}`$3$`d^9`$ components, but only the 2$`p_{3/2}`$3$`d^{10}`$ peaks are shown in the figure. One can see the obvious doping dependent shifts of O $`1s`$ and Nd $`3d`$ core levels from both the displaced and overlayerd plots in Fig. 1. To deduce the amount of the core-level shifts reliably, we used the peak position for the Nd 3d spectra and the mid point of the lower binding energy slope for the O 1s spectra. We used the mid-point position rather than the peak position for O $`1s`$ because the line shape on the higher binding energy side of the O 1s peak was sensitive to a slight surface degradation or contamination. The Cu 2p core-level line shape was not identical between different $`x`$’s, and becomes broader as $`x`$ increases. This is because the doped electrons in the CuO2 plane produce Cu1+ sites on the otherwise Cu2+ background, which yield an overlapping chemically shifted component located on the lower binding side of the Cu2+ peak. Therefore, it was difficult to uniquely determine the shift of the Cu 2p core level and we only take its peak positions in the following.

Figure 2 shows the binding energy shift of each core level relative to the as-grown $`x=0.05`$ sample. Here, we have assumed that the change of the electron concentration caused by the oxygen reduction was $`0.04`$ per Cu (oxygen reduction being $`0.02`$) as reported previously . One can see that the Nd 3d and O 1s levels move toward higher binding energies with electron doping. The shift of Cu 2p is defined by the shift of the peak position, and is in the opposite direction to Nd 3d and O 1s because of the Cu1+ components mentioned above. We also measured the shifts of the core levels using the Mg $`K\alpha `$ line and almost the same results were obtained as shown in Fig. 2.

While the shift of the chemical potential changes the core-level binding energy, there is another factor that could affect the binding energy, that is, the change in the Madelung potential due to Ce4+ substitution for Nd3+. However, the identical shifts of the O 1s and Nd 3d core levels indicate that the change in the Madelung potential has negligible affects on the core-level shifts because it would shift the core levels of the O2- anion and Nd3+ cation in the opposite directions. Moreover, as the shifts of the O $`1s`$ and Nd $`3d`$ core levels toward higher binding energies with electron doping are opposite to what would be expcted from increasing core-hole screening capability with $`x`$, excluding the core-hole screening mechanism as the main cause of the core-level shifts. Therefore, we conclude that the shifts of the O $`1s`$ and Nd $`3d`$ core levels are largely due to the chemical potential shift $`\mathrm{\Delta }\mu `$. We have evaluated $`\mathrm{\Delta }\mu `$ in NCCO by taking the average of the shifts of the two core levels.

Figure 3(a) shows $`\mathrm{\Delta }\mu `$ in NCCO as well as $`\mathrm{\Delta }\mu `$ in LSCO  as a function of electron or hole carrier concentration. In order to obtain the jump in $`\mu `$ between Nd2CuO4 and La2CuO4, we also measured the O $`1s`$ and Cu $`2p`$ levels in LSCO as shown in Fig. 4, and found that the O $`1s`$ and Cu $`2p`$ levels in Nd2CuO4 lie at $`150`$ meV and $`400`$ meV higher binding energies than those in La2CuO4, respectively. The fact that the observed jump is different between O 1s and Cu 2p is not surprising because Nd2CuO4 and La2CuO4 are different materials with different crystal structures. Thus the chemical potential jump between Nd2CuO4 and La2CuO4 cannot be uniquely determined from those data but it should be much smaller than the CT gap of about $`1.5`$ eV for Nd2CuO4 and $`2.0`$ eV for La2CuO4 estimated from the optical measurements . This small jump is in accordance with the early valence-band photoemission studies of LSCO and NCCO .

Figure 3(a) demonstrates the different behaviors of $`\mathrm{\Delta }\mu `$ between LSCO and NCCO. In LSCO, $`\mathrm{\Delta }\mu `$ is suppressed in the underdoped region $`x0.13`$, whereas in NCCO $`\mathrm{\Delta }\mu `$ monotonously increases with electron doping in the whole concentration range. Figure 3(c) represents the phase diagram of LSCO and NCCO drawn against the chemical potential $`\mu `$ and the temperature $`T`$. One can see that the $`\mu T`$ phase diagram is rather symmetric between the hole doping and electron doping unlike the widely used $`xT`$ phase diagram. That is, in both the electron- and hole-doped cases, the SC region is adjecent to the AF region, as proposed by Zhang  based on SO(5) symmetry. The present phase diagram implies that as a function of $`\mu `$, $`T_c`$ increases up to the point where the superconductivity is taken over by the AF ordering. Such a behavior is reminiscent of the superfluid-solid transition in the $`pT`$ phase diagram of 3He .

The monotonous increase of $`\mu `$ in NCCO may be understood within the simple rigid-band model. The chemical potential of Nd2CuO4 lies near the bottom of the conduction band. When electrons are doped, $`\mu `$ moves upward into this band as long as the AF ordering and hence the AF band structure persist as in NCCO. This behavior is contrasted with the suppression of $`\mathrm{\Delta }\mu `$ in underdoped LSCO, where the AF ordering is quickly destroyed and the electronic structure is dramatically reorganized by a small amount of hole doping. It has been suggested  that the suppression of $`\mathrm{\Delta }\mu `$ is related to the charge fluctuations in the form of stripes in LSCO. Indeed in La2-xSrxNiO4 (LSNO), where static stripe order is stable at $`x`$ 0.33 , $`\mathrm{\Delta }\mu `$ is anomalously suppressed in the underdoped region $`x0.33`$ . In LSCO, dynamical stripe fluctuations have been implied by inelastic neutron scattering studies  whereas any sign of stripes is absent in the neutron study of NCCO . In Fig. 3, we compare the chemical potential shift $`\mathrm{\Delta }\mu `$ with the incommensurability $`ϵ`$ which was deduced from the neutron experiments both for LSCO  and NCCO . This figure shows that in the region where $`\mu `$ does not move ($`x0.13`$ in LSCO), $`ϵ`$ linearly increases with $`x`$. In NCCO, where the chemical potential $`\mu `$ monotonously moves upward, the incommensurability does not change with $`x`$ (remains zero), in other words, static nor dynamical stripes do not exist. In the overdoped region of LSCO, the number of stripes saturates, doped holes overflow into the interstripe region and $`\mu `$ moves fast with hole doping.

The smallness of the chemical potential jump between Nd2CuO4 and La2CuO4 indicates that $`\mu `$ lies within the CT gap of Nd2CuO4 ($``$1.5 eV) and La2CuO4 ($``$2.0 eV). This behavior was clearly observed in an angle-resolved photoemission (ARPES) study of LSCO , where $`\mu `$ is located $``$0.4 eV above the top of the valence band of La2CuO4. Such a behavior is quite peculiar from the view point of the rigid-band model, and therefore a dramatic change in the electronic structure should occur with hole doping. In NCCO, however, the doping-induced change is not so dramatic as in LSCO as mentioned above, and the chemical potential pinning well below the bottom of the conduction band is not very likely. Another possible cause of the small chemical potential jump is that the CT gap is indirect and is smaller than that estimated from the optical studies. In optical conductivity spectra, only the direct transition can be unambiguously measured, and if the gap is indirect, the gap would be estimated much larger than the CT gap. This idea is consistent with the recent resonant inelastic x-ray scattering study  and its theoretical analysis using $`t`$-$`t^{^{}}`$-$`t^{^{\prime \prime }}`$-$`U`$ model , which has yielded an indirect gap that is smaller than the direct one by $`0.5`$ eV.

In summary, we have experimentally determined the doping dependence of the chemical potential shift $`\mathrm{\Delta }\mu `$ in NCCO and observed a monotonous shift with doping. Comparison with LSCO indicates that the change in the electronic structure with carrier doping is more moderate in NCCO. The monotonous shift is consistent with the observation that spin fluctuations are commensurate in NCCO. The small chemical potential jump between the $`n`$-type and $`p`$-type materials is confirmed and is attributed to the indirect CT gap and the chemical potential pinning within the CT gap in LSCO.

The authors would like to thank S. Tesanovic, A. Ino and T. Mizokawa for enlightening discussions. Collaboration with G. A. Sawatzky and J. van Elp in the early stage of this work is gratefully acknowledged. This work was supported by a Grant-in-Aid for Scientific Research in Priority Area “Novel Quantum Phenomena in Transition Metal Oxides” and a Special Coordination Fund for the Promotion of Science and Technology from the Ministry of Education and Science and by New Energy and Industrial Technology Development Organization (NEDO).


format: text

text:
# The Promise of String Cosmology 11footnote 1This article was commissioned in conjunction of the Workshop on String Cosmology held from July 24 - Aug. 4 2000 at PIMS, Univ. of British Columbia and is published in APCTP Bulletin # 6, (2001). Brown preprint BROWN-HET-1243.

## Abstract

The interplay between string theory and cosmology is very promising. Since string theory will yield a quantum theory of space-time and unify all forces of nature, it has the potential of addressing many of the conceptual problems of the current models of early Universe cosmology. In turn, cosmology is the most obvious testing ground in the effort to construct non-perturbative string theory, and can provide the crucial connection between theory and experiment/observation 222For a more technical overview of the existing approaches to connecting string theory and cosmology, the reader is referred in ..

Early Universe cosmology has made spectacular progress over the past two decades, driven both by a wealth of new observational data and by theoretical breakthroughs. In particular, a scenario for the very early Universe, the inflationary Universe scenario , has emerged which has led to solutions of some of the deep mysteries of standard big bang cosmology, and which has made some quite specific predictions for the geometry of the Universe and for the spectrum of cosmic microwave background (CMB) anisotropies which have been confirmed by the most recent observations .

In spite of its success, the inflationary Universe scenario suffers from some deep conceptual problems  which have prevented it from achieving the status of an honest physical theory. For example, scalar-field-driven inflationary cosmology is incomplete in the same sense as the standard big-bang cosmology was: in the context of inflationary cosmology an initial singularity is inevitable, thus making it impossible to formulate a consistent initial value problem. More seriously, the successful predictions of inflationary cosmology for the spectrum of the CMB depend - in the current models of inflation in which the exponential expansion of space is driven by the potential energy of a scalar matter field - on extrapolating the physical theory into a regime where its foundations (namely classical general relativity and weakly coupled quantum field theory) are known to break down 333This is the so-called trans-Planckian problem for inflationary cosmology .. The Achilles heel of inflationary cosmology is the cosmological constant problem. How do we know that the unknown mechanism which cancels the vacuum energy today does not also cancel the transient vacuum energy which is required to drive inflation? Quite apart from these conceptual problems, a convincing realization of inflation in the context of our current models of particle physics based on quantum field theory is lacking.

Recent observations are also beginning to challenge some of the premises of our current cosmological models. These observations point to the existence of a dark energy component in the Universe which makes up the bulk of the energy density but which does not cluster . Furthermore, recent measurements and numerical simulations show a mismatch between the predicted and observed galaxy halo structures which may imply that the dark matter is not cold (see e.g.  and references therein).

String (or M-) theory is the best candidate for a unified theory of all forces and as such for a quantum theory of space-time and matter. Thus, it will resolve the conflict between quantum field theory and general relativity, and therefore it has the potential of addressing many of the open issues in early Universe cosmology. If successful in this challenge, string theory will have made an important contact with experiment/observation, and will no longer be able to be criticized for lacking a connection to data and for thus being pure mathematics rather than a physical science. The successful determination of the entropy of an extremal black hole in the context of string theory  is a first major success in establishing a link between string theory and cosmology.

Until recently there has not been much work on the interface between string theory and cosmology. However, it has been clear since the late 1980’s that string theory has the potential to solve some of the problems of standard cosmology. For example, it was argued , in the context of perturbative string theory, and assuming that all spatial dimensions are toroidally compactified, that as a consequence of t-duality the physics at very small radii is equivalent to that at very large radii, hence eliminating the cosmological singularities. Moreover, string winding modes may yield a mechanism which allows at most three spatial dimensions to become large, thus providing a dynamical resolution of the potential embarrassment that superstring theory is consistent only in nine (or ten in the context of M-theory) spatial dimensions. This example demonstrates that with string theory, questions about the Universe can now be addressed scientifically which before the advent of string theory were exclusively in the realm of philosophy.

In the past five years string theory has undergone a “second revolution”. With the discovery of D-branes  it has become clear that there are many more fundamental degrees of freedom than just strings. Dualities between the previously known five consistent string theories were discovered which indicate that there is a common underlying theory, M-theory. The different string theories (and also 11-d supergravity) correspond to certain corners of moduli space of M-theory. A further new important tool is the AdS-CFT correspondence  which relates classical gravity on an anti-de-Sitter background to a conformal field theory living on the boundary. Recently, much attention has been focused on brane world scenarios, scenarios based on a combination of the above ideas in which the matter fields of our Universe are confined to a three-brane and only gravity lives in the higher-dimensional bulk (see e.g. ). In general, the spatial dimensions transverse to the branes must be very small in order to avoid too large deviations from Newton’s gravitational laws. An interesting suggestion first put forward without a direct connection to string theory  444See  for a realization of this idea in the context of string theory. is that two of the extra dimensions are large (mm scale). Since the four dimensional effective gravitational constant is related to the fundamental higher dimensional one via the size of the extra dimensions, the fundamental gravitational constant can be reduced to the TeV scale, provided the large extra dimensions have a size of the order of 1mm. This provides the potential for solving the mass hierarchy problem of elementary particle physics. By constructing models in which gravity remains confined to the brane, it is possible to make the extra dimensions large (in fact even infinite ).

On the other hand, there is still no non-perturbative formulation of M-theory; not even its fundamental physical degrees of freedom are known. In at least one approach to M-theory , the underlying theory is postulated to be a quantum mechanical matrix model. If this is true, then some of the non-commutativity of the matrix model will be reflected in non-commutativity of the space-time coordinates at a microscopic level. Returning to the more traditional description of M-theory, another serious problem results from the fact that the moduli space of vacuum states corresponding to any of the string theories is huge, and many of the predictions for particle physics as well as for cosmology depend on which vacuum state is chosen. A key challenge to string cosmology is to find a mechanism which distinguishes the vacuum state corresponding to our Universe. Since the moduli correspond to massless fields, new problems for cosmology appear: why don’t we see the massless particles associated with the moduli fields, and if this question is answered in the usual way by invoking some non-perturbative mass generation for the moduli, then why don’t the moduli over-close the Universe ?

However, even without first having to solve these difficult problems, there are interesting avenues for string cosmology. A particularly interesting issue is to determine whether there are string-specific ways to obtain a period of inflationary expansion in the early Universe. If the answer to this question is affirmative, then it is of utmost importance to try to find observational signatures which could distinguish these string-specific inflationary models from scalar-field-driven inflation. The answer to this question will inevitably require looking beyond the power spectra of density fluctuations and of CMB anisotropies, since any set of power spectra can be obtained by suitable choice of a toy model scalar field potential. It may, however, also turn out that string cosmology does not yield anything that looks like the usual inflationary scenario, but that the resulting theory nevertheless provides an alternative to inflation in solving the problems of standard big bang cosmology such as the horizon and flatness problems. This possibility is realized in the pre-big-bang approach to string cosmology  based on dilaton gravity as the low-energy effective theory consistent with the symmetries of string theory. In this scenario the background in the Einstein frame is in fact contracting super-exponentially in the pre-big-bang phase, leaving an effective Hubble constant whose absolute value is increasing in time. In pre-big-bang cosmology the spectrum of scalar metric fluctuations stemming from the dilaton gravity sector is not scale-invariant as it would be in inflationary cosmology . A more radical alternative to inflation which can solve the homogeneity and flatness problems is the varying-speed-of-light scenario in which the speed of light in the early stages of the evolution of the Universe is postulated to be many orders of magnitude larger than at the present time .

Another direct observational issue which string cosmology should address concerns the origin of dark energy. Since string theory leads to a large number of scalar moduli fields, it is of great interest to explore the potential of these fields for explaining the origin of dark energy. Possibly related to this issue is the question of the smallness of the observed cosmological constant. Finally, even without knowing the exact vacuum state of string theory it should be possible to study whether string theory will indeed wash out the singularities of standard and inflationary cosmology.

Interest in string cosmology is growing. From July 24 - August 4, 2000, about 45 physicist working on string theory, non-perturbative field theory and cosmology gathered at the University of British Columbia (UBC) in Vancouver for a workshop devoted to string cosmology. This workshop was co-sponsored by the Pacific Institute for Mathematical Sciences (PIMS), the Canadian Institute for Advanced Research (CIAR) and by the APCTP. The invited participants included B. Greene (Columbia Univ.), N. Kaloper (Stanford), L. Kofman (CITA, Univ. of Toronto), D. Lowe (Brown Univ.), B. Ovrut (Univ. of Pennsylvania), S. Ramgoolam (Brown Univ.), S. Sin (Hanyang Univ.), D. Son (Columbia Univ.), P. Steinhardt (Princeton Univ.), H. Verlinde (Princeton Univ.), G. Veneziano (CERN) and A. Zhitnitsky (UBC) 555See the web site http://kepler.physics.ubc.ca/p̃fs99 for details about the workshop.. Many different approaches to string cosmology were represented at this workshop. In the following discussion of various approaches to string cosmology I will limit myself to those which were discussed extensively at the workshop.

Based on the conclusions of one of the most conservative approach to string cosmology, namely Pre-Big-Bang Cosmology , it is clear that we should expect important differences in the evolution of the early Universe in string cosmology compared to scalar-field-driven inflationary cosmology. Pre-big-bang cosmology is a low energy description of physics which takes into account all of the low energy modes in string theory, i.e. not only the graviton but also the dilaton and the anti-symmetric tensor field (the latter, however, does not play an important role for the background dynamics at the level of the homogeneous equations of motion 666As was mentioned at the workshop, a non-vanishing tensor field decreases the set of inhomogeneous initial conditions which can lead to the pre-big-bang evolution.). The symmetries and equations of motion for dilaton gravity give rise to a scenario of cosmology in which the Universe starts out in a state near the perturbative string vacuum, goes through a dilaton-driven phase of contraction during which the Hubble radius $`H^1(t)`$ ($`H`$ being the expansion rate of the Universe and $`t`$ denoting time) shrinks faster than the physical length of fixed comoving scales, thus simultaneously solving the horizon problem of standard cosmology and giving rise to a mechanism for structure formation similar as in the usual inflationary cosmologies: quantum vacuum fluctuations existing during the contraction phase on sub-Hubble lengths get frozen in when the Hubble radius becomes smaller than the wavelength of the fluctuation, they are amplified super-adiabatically while the wavelength is larger than the Hubble radius, and then re-enter the Hubble radius in the late Universe as classical density fluctuations. As the Universe contracts, the curvature increases. Eventually, in a way which unfortunately cannot be described by the equations of dilaton gravity, the Universe makes a transition to an expanding Friedmann cosmology related to the contracting phase via a duality transformation. Note that in pre-big-bang cosmology, both the background evolution and the spectrum of induced fluctuations is different from that of scalar-field-driven inflationary models. The density fluctuations from the dilaton sector are not scale invariant. In fact, they are completely unimportant on large length scales. The other fields present in pre-big-bang cosmology (such as the axion) will, however, generate fluctuations which can be scale-invariant . Since these fluctuations start out as isocurvature fluctuations, they lead to different predictions for observables such as the spectrum of CMB anisotropies than the scalar metric fluctuations of inflationary cosmology. The lesson we learn is that it is quite likely that the predictions of string cosmology will differ enough from those of standard inflationary cosmology to enable observations to probe the physics of the Planck scale.

Pre-big-bang cosmology suffers from the graceful exit problem: in the context of the equations of motion of dilaton gravity, the contracting branch of the evolution terminates in a singularity, and the expanding branch emerges from an equivalent singularity. Since the singularity occurs in the high curvature regime, it is clear that the effective dilaton gravity action no longer describes the stringy physics in a correct way close to the singularity. A challenge for string cosmology is to understand if true stringy physics can solve the graceful exit problem and lead to a smooth transition from the contracting dilaton-dominated phase to the expanding phase with fixed value of the dilaton. Although there have been several attempts to go beyond dilaton gravity in order to solve the graceful exit problem (see e.g.  for a recent overview), there is not yet a convincing approach. An important open problem for string cosmology is to develop a truly stringy, non-perturbative version of the pre-big-bang scenario, or to find out how the scenario must be modified in order to become a truly stringy cosmological model.

Another promising approach to string cosmology is based on Horava-Witten theory  and is sometimes called heterotic M-theory cosmology . The theory is based on one particular ray in the moduli space of M-theory, namely the compactification of 11-dimensional supergravity on an orbifold $`S^1/Z_2`$. This leads to two distinguished co-dimension one planes, the orbifold fixed planes. As usual, it is assumed that, in addition, six spatial dimensions are compactified on a special type of Calabi-Yau manifold, and that there is a set of $`E_8`$ gauge supermultiplets on each orbifold plane. The reason for adopting this particular orbifold compactification is that it allows for chiral fermions . The resulting theory is dual to the strongly coupled $`E_8\times E_8`$ heterotic string theory. Matching at tree level to the phenomenological gravitational and grand unified gauge couplings, one finds that the orbifold must be larger than the Calabi-Yau radius. This suggests that there is a regime in which the dynamics of the Universe is effectively five-dimensional. The resulting five-dimensional effective theory is a gauged $`N=1`$ supergravity coupled to four-dimensional boundary gauge and matter supermultiplets (see e.g.  for a recent survey). Initially, only the effective four-dimensional equations were studied. However, it has recently been realized that the dynamics in the orbifold direction is quite nontrivial and may lead to new effects of great interest to cosmology. In particular, the equations of motion in five dimensions admit as a static solution a pair of three-branes parallel to the orbifold fixed planes. It is tempting to try to obtain inflation from the nontrivial dynamics involving the orbifold direction. Note, however, that as in pre-big-bang cosmology, many truly stringy effects are frozen out in the effective action-based analysis of heterotic M-theory cosmology.

A rather different approach to string cosmology is based on proceeding in analogy to the usual starting point of big bang cosmology, which is the assumption that the Universe starts out with Planck size and emerges from a hot initial state in which all degrees of freedom were highly excited. It seems quite natural to generalize this starting point to string cosmology and to assume  that the Universe starts out with all spatial dimensions of string scale and all fundamental degrees of freedom excited according to their thermal distribution. In particular, this means that the winding modes of strings and branes are highly excited. This approach is called the brane gas approach to string cosmology. In this approach, all spatial dimensions start out compact and the question for cosmology is not why certain dimensions are compactified, but why a fixed number of dimensions can dynamically decompactify (in the sense of their radii becoming larger than the present Hubble radius of the Universe). There have been some interesting recent developments  along the lines of this approach to string cosmology. In particular, assuming that all spatial dimensions are compactified on a torus, it has been shown that whereas at early times (small sizes) the higher dimensional branes dominate the evolution, it is the fundamental strings which have the most important effect at late times . Thus, the mechanism of , which only allows three spatial dimensions to become large, survives in the modern view of string theory. In addition, the duality under inversion of the radius of the torus is maintained, thus assuring that no physical cosmological singularities will arise. It is important to generalize these considerations to more realistic background spaces.

The approach to string cosmology which at the present time is attracting most attention is the class of brane-world scenarios (see e.g.  for work along these lines presented at the APCTP/CIAR/PIMS string cosmology workshop). The starting point of this approach is the assumption that our matter fields are confined to a three-brane (four space-time dimensions) embedded in a higher-dimensional (usually five-dimensional) bulk space-time. An important realization is that it is possible  to confine gravity on one of the branes by having the four-dimensional metric depend on the fifth coordinate (the orbifold coordinate in the language of heterotic M-theory cosmology). In many cases, the bulk metric is AdS, and this gives rise to a tantalizing connection with the AdS-CFT correspondence. However, most work along these lines is not directly justified from string theory. There appear to be an unlimited number of versions of the basic scenario. As was stressed at the APCTP/CIAR/PIMS string cosmology workshop, an urgent task for string cosmology is to derive the basic rules which a brane-world scenario must obey if it is to emerge from string theory. Is the brane an orbifold fixed plane (as it is in heterotic M-theory cosmology) or is it a D-brane? How many branes are there in the scenario? Are there matter fields in the bulk? At the present time, most work on brane-world cosmology is based on studying the background gravitational equations of motion of the brane in the presence of the bulk, making use of the general relativistic Israel matching conditions. It has been realized that the equations for the dynamics of our brane in the brane-world scenario may be quite different from the usual four-dimensional general relativistic equations. This provides stringent constraints on brane world cosmologies. However, it also provides the chance to obtain new effects, or to provide new sources of inflation. In order to make contact with the wealth of recent cosmological observations on the distribution of light and CMB light, it is crucial to extend the theory of cosmological perturbations, the basis for making such a contact, from the usual four-dimensional setting to the brane-world scenarios. Initial work along these lines was presented at the workshop (see  and references therein). A limitation which the brane-world scenarios share with some of the other approaches to string cosmology is that most of the analysis is based on a truncation of the dynamics to low energy modes, thus eliminating many truly stringy effects. It is, however, possible to obtain interesting new effects for cosmology. For example, if our Universe is a three-brane moving in a bulk black hole background, it is possible  to find a realization of the varying light scenario for solving the problems of standard cosmology.

Ultimately, however, a complete cosmological model should be based on the non-perturbative formulation of M-theory. Starting from the matrix theory approach to M-theory (see e.g. ), it seems likely that the resulting space-time theory will be non-commutative at the string scale (see e.g.  for work along these lines presented at the workshop). One limitation of the matrix theory approach to M-theory  is that the theory is at the present time only defined in light-cone coordinates, which makes the connection to a space-time formulation very challenging.

A generic problem which appears in all approaches to string cosmology is the moduli problem. What freezes the value of the dilaton and what sets the scale of the extra dimensions? Usually it is assumed that some non-perturbative physics is responsible for the answers. If an effective potential provides the required mechanism (see e.g.  for some recent work presented in Vancouver), then why does the value of the potential at its minimum not lead to another source for the cosmological constant? Thus, the moduli problem and the cosmological constant problem appear to be related. There are some new approaches to solving the cosmological constant problem based on the brane-world scenario, e.g. a model by Verlinde  in the context of a theory with two branes, the physical brane and the Planck brane, which assumes that the cosmological constant on the Planck brane is protected by supersymmetry, and that it flows to a small cosmological constant on the physical brane via renormalization group flow.

The field of string cosmology is still in its infancy. Although many interesting approaches to connecting string theory and cosmology exist, most of them suffer from our ignorance about the true non-perturbative formulation of string theory. The promises of the field, however, are great, as witnessed by the excitement expressed at the APCTP/CIAR/PIMS string cosmology workshop

Acknowledgements

I wish to thank the APCTP, CIAR and PIMS for their generous support of the String Cosmology Workshop held at UBC in Vancouver as part of the Frontiers of Mathematical Physics workshop series jointly organized by PIMS and APCTP. I also thank R. MacKenzie for useful comments on the draft. The author’s research was supported in part by the U.S. Department of Energy under Contract DE-FG02-91ER40688, TASK A.


text:
# Light Curves from an Expanding Relativistic Jet

## 1 The Hydrodynamics

The hydrodynamic simulations are performed using a fully relativistic 2D AMR code. The jet propagates into a cold and homogeneous ambient medium of number density $`n_1\mathrm{cm}^3`$. We use $`n_1=1`$, typical of the ISM. We assume an adiabatic evolution (i.e. that radiation losses do not influence the dynamics). For the initial conditions, we use a wedge of opening angle $`\theta _0=0.2`$ taken out of the Blandford McKee BM  self similar spherical solution, with an isotropic equivalent energy of $`E_{52}10^{52}`$ ergs, where we use $`E_{52}=1`$. Since the lateral expansion is expected to become significant only when the Lorentz factor of the flow drops to $`\gamma 1/\theta _0=5`$, the initial Lorentz factor of the shock was chosen to be $`\mathrm{\Gamma }=\sqrt{2}\gamma =23.7`$.

Figure (1) shows a 3D view of the jet at the last time step of the simulation, with color-maps of the proper number density, $`n^{}`$, and proper emissivity. Primed quantities are measured in the local rest frame of the fluid. While the number density does not change significantly between the front and the sides of the jet, the emissivity is significant only at the front of the jet, and drops sharply at angles larger than the initial opening angle, $`\theta _0`$. The overall egg-shaped structure is very different from the quasi-spherical structure assumed in 1D analytic models.

## 2 The Light Curves

The synchrotron spectral emissivity, $`P_\nu ^{}^{}`$, is calculated using the physical conditions determined by the hydrodynamic simulation. The electrons and the magnetic field are assumed to hold fractions $`ϵ_e`$ and $`ϵ_B`$, respectively, of the internal energy density, $`e^{}`$, while the electrons posses a power law energy distribution, $`N(\gamma _e)\gamma _e^p`$ for $`\gamma _e>\gamma _m=\left[(p2)/(p1)\right](ϵ_ee^{}/n^{}m_ec^2)`$. For simplicity, we ignore the effects of electron cooling and self absorption, so that $`P_\nu ^{}^{}\nu ^{1/3}`$ for $`\nu ^{}<\nu _m^{}`$ and $`P_\nu ^{}^{}\nu ^{(1p)/2}`$ for $`\nu ^{}>\nu _m^{}`$, where $`\nu _m^{}`$ is the local synchrotron frequency of an electron with $`\gamma _e=\gamma _m`$; $`F_\nu `$ is calculated using the formalism of GPS , summing over the contributions from the finite 4-volume of the simulation.

Figure (2) shows the radio light curves seen by an observer along the jet axis ($`\theta _{obs}=0`$). For simplicity, cosmological corrections are not included. The insert shows an optical light curve as seen by observers at three different viewing angles with respect to the jet axis: $`\theta _{obs}/\theta _0=0,1,2`$. We obtain an achromatic ‘jet break’ in the light curve for $`\nu >\nu _m(t_{jet})`$ (as predicted by simple semi-analytic models Rhoads ; SPH ) at $`t_{jet}(\theta _{obs})`$, where $`0.66t_{jet}(\theta _0)t_{jet}(0)5.8(E_{52}/n_1)^{1/3}(\theta _0/0.2)^{8/3}`$ days. Defining $`\alpha `$, $`\beta `$ by $`F_\nu \nu ^\beta t^\alpha `$, the shape of the break may be approximated by

$$F_\nu =F_0\nu ^\beta \left[(t/t_{jet})^{s\alpha _1}+(t/t_{jet})^{s\alpha _2}\right]^{1/s},$$
(1)
where $`\alpha _1=\alpha (tt_{jet})`$, $`\alpha _2=\alpha (tt_{jet})`$ and the parameter $`s(\theta _{obs})`$ determines the sharpness of the break, and ranges between $`s(\theta _0)1`$ to $`s(0)4.5`$, indicating that the break is sharper at smaller $`\theta _{obs}`$. For $`\nu <\nu _m(t_{jet})`$, there is only a moderate and more gradual change in $`\alpha `$, until the time when $`\nu _m=\nu `$. For $`\nu >\nu _m`$ we find that $`\alpha _2`$ is slightly smaller than $`p`$ (the value predicted by most simple models) for $`\theta _{obs}=0`$ ($`\alpha _2=2.85`$ for $`p=2.5`$).

## 3 Discussion

We find that the physical conditions at the sides of the jet are significantly different than at the front of the jet, and most of the radiation is emitted within the initial opening angle of the jet \[$`\theta <\theta _0`$; see Figure (1)\]. Therefore, the frequently used assumption of a homogeneous jet seems inadequate. For $`\nu >\nu _m(t_{jet})`$ we find a sharp achromatic break in the light curve at $`t_{jet}=5.8(E_{52}/n_1)^{1/3}(\theta _0/0.2)^{8/3}`$ days, for $`\theta _{obs}=0`$, while at larger viewing angles $`t_{jet}`$ increases and the break becomes smoother. For $`\nu <\nu _m(t_{jet})`$ the change in the temporal index $`\alpha `$ near $`t_{jet}`$ is more moderate and gradual. The value of $`\alpha `$ for $`t>t_{jet}`$, $`\nu >\nu _m`$, $`\theta _{obs}=0`$ is slightly smaller than $`p`$ ($`\alpha =2.85`$ for $`p=2.5`$). Finally, we note that $`t_{jet}(\theta _{obs}=\theta _0)1.5t_{jet}(\theta _{obs}=0)`$. Since, in order to detect a burst in $`\gamma `$-rays we require that $`\theta _{obs}<\theta _0`$, this may induce an uncertainty of up to $`15\%`$ when deducing the value of $`\theta _0`$ from $`t_{jet}`$, unless $`\theta _{obs}`$ is well constrained.


text:
# Acknowledgement

## Acknowledgement

I wish to thank Michael Müller-Preussker and his group for their kind hospitality.


text:
# Magnetic Stresses at the Inner Edges of Accretion Disks Around Black Holes

## I Introduction

Ever since the seminal work of Novikov and Thorne jhk:nt73  and Shakura and Sunyaev jhk:ss73 , our understanding of the physics of accretion disks around black holes has been built on two fundamental beliefs: the energy available for radiation is equal to the binding energy of the innermost stable orbit, and is therefore determined solely by the spin of the black hole; and the (vertically-integrated) stress responsible for transporting angular momentum outwards through the disk is proportional to the local pressure. Recent work involving magnetic forces in the innermost regions of relativistic accretion disks has undercut both of these claims. In this review I will summarize the reasons for this change of heart and explain the consequences that are likely to follow if the new point of view proves correct.

## II The Traditional View of the Inner Edges of Relativistic Accretion Disks

The basic framework for understanding the inner regions of relativistic accretion disks was laid out more than twenty-five years ago jhk:nt73 ; jhk:pt74 . In the simplest possible picture (i.e., one assuming that disks are time-steady, axi-symmetric, and geometrically thin), their radial structure may be defined in terms of conservation equations that appear almost Newtonian, as the general relativistic effects can be collapsed into multiplicative correction factors. Conservation of angular momentum, for example, may be expressed through a single first-order differential equation in radius.

Thus, to specify the entire solution, all that is necessary is to choose a boundary condition for this differential equation. This boundary condition may be physically interpreted as determining the conserved outward angular momentum flux through the disk, or, equivalently, as specifying $`T_{r\varphi }`$, the $`r`$$`\varphi `$ component of the stress exerted on the disk, at its inner edge. This inner edge is conventionally taken to lie at the radius of the marginally stable orbit, $`r_{ms}`$.

Selecting this boundary condition is the only part of the procedure in which there is any guesswork. Until very recently, it was almost universally assumed that the stress must be zero at $`r_{ms}`$. Two reasons were commonly cited for this choice. The first, due to Thorne and collaborators jhk:nt73 ; jhk:pt74 , was that the radial speed inside $`r_{ms}`$ is so much larger than the radial speed outside $`r_{ms}`$ that, assuming a constant mass accretion rate, the inertia of mass inside $`r_{ms}`$ must be much smaller than the inertia of the disk proper in the region of stable orbits. Given that, it appears difficult to see how the small amount of matter in the plunging region could exert any significant force on the much heavier disk. The second argument jhk:ak89  begins with the assumption that the stress carrying angular momentum outward should always be proportional to the local pressure; i.e., $`T_{r\varphi }=\alpha p`$jhk:ss73 . If this is so, as material accelerates inward and expands, its pressure falls rapidly, so the stress must do likewise.

If there is no stress inside $`r_{ms}`$, then no forces change the energy or angular momentum of material inside that point. Matter therefore arrives at the event horizon with the energy and angular momentum it had at $`r_{ms}`$. The maximum energy per unit mass available to be radiated in the disk (the maximum radiative efficiency) is then defined simply by the binding energy of the marginally stable orbit, and the spin-up rate of the black hole is the mass accretion rate times the specific angular momentum of that orbit. Both the specific binding energy and the specific angular momentum at $`r_{ms}`$ depend only on the black hole’s normalized angular momentum $`a/M`$; thus, it has long been thought that both the radiative efficiency and the spin-up rate per unit accreted mass are functions only of $`a/M`$.

Interestingly, the very first work on this subject jhk:pt74  acknowledged that strong magnetic fields could upset the arguments that the stress must go to zero at $`r_{ms}`$, but chose (reasonably, given the state of knowledge at the time) to ignore this possibility. However, work of the past decade (summarized in jhk:bh98 ) has shown that magnetic fields are essential to angular momentum transfer in accretion disks. We must therefore re-open the question of the appropriate boundary condition on $`T_{r\varphi }`$ at $`r_{ms}`$.

## III Magnetic Fields in Accretion Disks

This recent work has shown that, in the main body of the disk, MHD turbulence creates a stress that, when vertically-averaged, is roughly proportional to the pressure, with an $`\alpha 0.01`$ – 0.1 jhk:bh98 ; jhk:s96 ; jhk:b96 . The magnetic stress accounts for the bulk of the angular momentum flux so long as two conditions are met: there must be some seed magnetic field in the accreting gas, but in the disk midplane the energy density of any imposed magnetic field must be smaller than the total pressure; and the matter must have high enough conductivity that the MHD approximation is valid. In the conditions surrounding black holes, these are all almost certain to apply. Weak magnetic fields are virtually ubiquitous in astrophysical environments; and near a black hole the gas temperature will very nearly always be high enough to maintain the gas in an ionized, and therefore highly conducting, state.

This turbulence is generated by an MHD instability driven by the orbital shear jhk:bh91 . Its growth rate is always comparable to the orbital frequency, so it grows roughly as rapidly (relative to an orbital period) in the relativistic part of the disk as in the Newtonian part. The nonlinear amplitude of the turbulence in the Newtonian part of the disk is determined by a cascade of energy toward shorter wavelengths, where it can ultimately be dissipated. This process, too, should operate more or less in the same fashion (as viewed in the fluid frame) in the relativistic part of the disk as in the non-relativistic part, provided the local inflow time is long enough for equilibrium to be established. We may therefore reasonably conclude that the ratio $`B^2/(8\pi p)`$ doesn’t change dramatically at radii a few times $`r_{ms}`$. Because orbital shear automatically stretches field lines in the sense that produces outward angular momentum flow ($`B_rB_\varphi <0`$), the effective $`\alpha `$-parameter ($`=B_rB_\varphi /(4\pi p)`$) is also unlikely to change much. But if this is so, why should the stress diminish as the marginally stable orbit is approached?

The only aspect of this process that changes in any qualitative way as $`r_{ms}`$ is approached from the outside is the ordering of four timescales: the inflow time $`t_{in}`$, the thermal time $`t_{th}`$, the turbulent dissipation time $`t_{diss}`$, and the dynamical time $`t_{dyn}`$. In the disk body, $`t_{in}\alpha ^1(r/h)^2t_{dyn}`$ and $`\beta t_{diss}t_{th}\alpha ^1t_{dyn}`$, where $`\beta `$ is the ratio of gas (+ radiation) pressure to magnetic energy density. Thus, $`t_{in}t_{th}>t_{diss}>t_{dyn}`$. However, just outside $`r_{ms}`$ the effective potential flattens (this is, of course, what it means for $`r_{ms}`$ to be the location of the innermost marginally stable orbit), so $`t_{in}`$ diminishes toward $`t_{dyn}`$. The concomitant decline in surface density likewise reduces $`t_{th}`$. However, $`t_{diss}`$ (crudely $`h/v_A`$ for disk thickness $`h`$ and Alfven speed $`v_A`$) declines more slowly than $`t_{th}`$. As a result, in the vicinity of $`r_{ms}`$, the interplay of plasma dynamics and flux-freezing should be at least as important to determining the field intensity as the balance between the turbulent dynamo and turbulent dissipation.

If flux-freezing really does determine the evolution of the field as the matter plunges inside $`r_{ms}`$, the magnetic field strength in the fluid frame stays roughly constant or increases somewhat even while the fluid pressure decreases dramatically jhk:k99 . It follows that the effective $`\alpha `$ rises sharply as the inflow accelerates near and inside $`r_{ms}`$. In fact, when the radial component of the velocity becomes relativistic, it is easy to show that these assumptions imply $`B^2/(8\pi )\rho c^2`$. That is, when the inflow is relativistic, magnetic forces become competitive with gravity, and the Alfven speed becomes relativistic. Matter may then retain a causal coupling with the disk it left behind even when it has fallen well inside $`r_{ms}`$, allowing significant transfer of energy and angular momentum jhk:k99 ; jhk:g99 .

Looking back at the old argument that the small inertia of mattter in the plunging region cannot do much to the main body of the disk, we now see that magnetic fields in effect turn this argument on its head. Magnetic connections between the disk and plunging matter act as “tethers” by which the massive disk restrains angular acceleration of the low-inertia streams inside $`r_{ms}`$. In so doing, stress can be exerted on the disk itself.

## IV Consequences

The existence of significant stress at the marginally stable orbit has major consequences for the most fundamental properties of accretion onto black holes. Energy taken from plunging matter and delivered to the disk is potentially available for radiation; removal of angular momentum from this matter retards the rate of black hole spin-up jhk:ak00 . In fact, when the black hole is rapidly rotating, the rotational energy of the black hole itself can be tapped: frame-dragging gives plunging matter a high orbital frequency; magnetic connections from this region to the disk exert a torque; the end-result is energy drawn from black hole rotation given to matter in the disk jhk:g99 ; jhk:k00 . In other words, the amount of energy drawn from accreting matter is no longer a quantity fixed by orbital mechanics: it is a dynamical quantity that is the product of complicated MHD dynamics, and may even be time-variable.

The fate of the work done by magnetic fields at $`rr_{ms}`$ depends critically on the ratio between the dissipation and inflow timescales in the place where the energy is brought. If the energy can be converted into heat and radiated before the matter finds its way into the plunging region inside $`r_{ms}`$, it adds to the radiative efficiency; if, on the other hand, dissipation is slow, so that the energy stays in the gas either as organized kinetic energy or as magnetic field energy, it may end up being taken into the black hole. Given the arguments of the previous section, it is as yet unclear how this balance works out.

If additional heat is deposited in the disk, there can be substantial alterations to the disk spectrum jhk:ak00 . Additional flux is radiated at high frequencies, and relativistic effects beam it strongly into the equatorial plane. In addition, because so much more radiation is produced in the most relativistic portion of the disk, the returning radiation fraction is greatly enhanced. This latter effect has implications for phenomenology as diverse as polarization of the emitted light and the synchronization of fluctuations.

However, there is an additional implication of stress at the inner edge of the disk that may be important even if the radiative efficiency is hardly altered by these mechanisms: dynamics at the marginally stable orbit can be a powerful “noise” source to the entire dynamical system. In the long run, these dynamical fluctuations may be reflected in luminosity fluctuations, which are, of course, a hallmark of all known accreting black holes. The origin of this “noise” may be seen from a simple thought experiment: Consider a small magnetized fluid element orbiting just outside $`r_{ms}`$. Imagine that half this fluid element receives a negative angular momentum perturbation from the turbulence and begins to fall inward. As it does so, its orbital frequency increases simply due to its decreasing orbital radius. The resulting shear between the two halves of the fluid element stretches the magnetic flux tube connecting them, and a magnetic tension force creates a torque that transfers angular momentum from the falling half to the half that is still orbiting stably. The stably orbiting half must then move outward, launching a compressive wave into the disk. Thus, the basic mechanics of accretion create disk “noise” when there are magnetic connections across the marginally stable orbit.

## V Simulations

All the arguments presented so far have been qualitative. Whether these effects are quantitatively important can only be ascertained as a result of genuine calculation. MHD turbulent dynamics are sufficiently complicated that (almost) all credible calculations are numerical simulations. Although much work remains to be done before simulations can be run with a requisite level of realism, some preliminary results have been achieved jhk:h00 ; jhk:hk01 ; jhk:arc01 .

To date, all simulations touching on these issues have been limited in a number of regards, both physical and numerical. All have assumed Newtonian physics; general relativistic dynamics is mimicked by the use of the Paczyński-Wiita potential $`U=GM/(r2GM/c^2)`$. This potential qualitatively reproduces motion in a Schwarzschild metric by creating a marginally stable orbit at $`r=6GM/c^2`$. In addition, all simulations so far have substituted an assumed equation of state for a real energy equation. This assumption has two deleterious consequences: magnetosonic waves don’t propagate at the right speed because the pressure isn’t correctly computed; and, more importantly, it is impossible to use these simulations to estimate the radiative efficiency because energy taken from the plunging region to the disk falls right back in. There are also aspects of the numerical methods used in these simulations that are less than optimal. Both kinetic energy and magnetic energy can disappear if fluid from two adjacent cells with oppositely directed velocity or magnetic field components is combined. It also appears likely that no simulation to date has used a grid fine enough to resolve the magnetic field structure.

Nonetheless, putting aside these qualms, the simulational results may still be used as an indication of what may appear in more realistic simulations. In particular, angular momentum transport should be more reliably treated than energy transport because it is not subject to the problems described in the previous paragraph. The result shown in Figure 1 jhk:hk01  is particularly noteworthy. Two distributions of $`T_{r\varphi }`$ integrated vertically and averaged over both azimuth and time are shown: the prediction of the Novikov-Thorne model for the mean mass accretion rate found in the simulation; and the magnetic contribution to the stress. The two curves coincide very closely over the radial range from $`8`$$`24GM/c^2`$, demonstrating that magnetic stresses account for very nearly all the angular momentum transport in the main part of the disk. They depart in an uninteresting way at large radius—the deviation here is due to the fact that the disk in the simulation had an outer edge, whereas the Novikov-Thorne model refers to a disk that extends to very large radius. At small radius, however, the contrast is very significant—as expected on the basis of the qualitative arguments presented earlier, magnetic stress remains important across the marginally stable orbit and throughout the region where gas plunges toward the event horizon. Indeed, as a result of this continuing stress, the mean specific angular momentum of matter crossing the inner edge of the simulation (at $`r=3GM/c^2`$) is about 5% smaller than the specific angular momentum at $`r_{ms}`$. Also as expected, the effective $`\alpha `$ parameter increases sharply in the innermost part of the accretion flow: between $`r=8GM/c^2`$ and $`r=4GM/c^2`$, it increases by an order of magnitude.

## VI Future Prospects

The results of these first simulations are encouraging, but they are far from the final answer. Fortunately, significant improvements are feasible. Improved resolution can be obtained both by cleverer gridding schemes and, of course, by a few years’ technological development. Real energy equations can be computed by explicitly incorporating phenemonological viscosity and resistivity and a radiative cooling function (following time-dependent radiation transfer is also possible jhk:ts01 , but will require somewhat greater computer power before it is feasible for this kind of simulation). Likewise, there is no fundamental impediment to writing MHD codes that work in a fixed relativistic metric. Thus, in a few years, we should be able to attach quantitative values to the effects pointed out qualitatively here.


text:
# A novel method of self interaction corrections in density functional calculations

## Introduction

Published as: Int. J. Quantum Chem. 81, 247 (2001).
Density functional theory (DFT) states that the ground state energy of a many particle system is an exact functional of the total particle density $`\rho `$. Unfortunately the exact functional is unknown and in order to have a practical mean of calculating the ground state energy one maps the interacting particle system onto a non-interacting, for which the total energy can be written as

$`E[\rho ]=T[\rho ]+E_\mathrm{H}[\rho ]+E_{\mathrm{XC}}[\rho ]+E_{\mathrm{NN}}+E_{\mathrm{eN}}[\rho ],`$ (1)
where $`T[\rho ]`$ is the one-electron kinetic energy (given by $`n_iϵ_iV_{\mathrm{eff}}\rho `$), $`E_\mathrm{H}[\rho ]`$ is the Hartree interaction. $`E_{\mathrm{XC}}[\rho ]`$, the exchange-correlation (XC) potential, is supposed to correct for the approximations introduced when replacing the true many-body kinetic energy by $`T[\rho ]`$ and the correct electron-electron interaction by $`E_\mathrm{H}[\rho ]`$. In addition $`E_{\mathrm{NN}}`$ is the electrostatic nucleus-nucleus interaction, and $`E_{\mathrm{eN}}[\rho ]`$ is the electron-nucleus interaction. Since the exact XC potential is not known, one has to rely to approximations. The simplest (and most frequently used) approximation is the LDA, which is a parameterized XC functional using a homogeneous electron density. Approximate expressions for the total energy combined with efficient methods for solving the Kohn-Sham equations have resulted in powerful methods for studying many materials properties such as magnetism, chemical bonding, hardness, elasticity as well as surface related phenomena, including relaxation, reconstruction and chemisorption. Quite frequently LDA calculations reproduce experimental results with an error smaller than a few percent. However, notable exceptions exist. One problem is that the approximations leading to the LDA introduce an unphysical self interaction (SI) of a particle with itself. For instance, in the local density approximation the Hartree part contains the mean-field interaction energy of an electron with itself. An exact XC functional would have a SI part which would exactly cancel the SI of the Hartree term. The result would be a theory totally free from self interaction, giving the exact ground state energy. It may be noted here that for systems with broad bands, where the electron states form Bloch states, self-interaction becomes neglibly small, whereas for systems where the electrons are localized or close to localization the self-interaction problem becomes serious. A self interaction correction (SIC) method of DFT was presented in 1981 by Perdew and Zunger (here referred to as PZ). Although there exist improved schemes the PZ-method is the method predominantly used. One of the virtues of the correction of PZ is that for the exact XC functional the added correction becomes zero. The method of PZ has been applied with some success to different systems, e.g. localized electrons in Ce and Pr, a 1-D Hubbard model and transition metal mono oxides. Although the PZ version of correcting for SI has been successful in many applications we question its ability to remove all self interaction and suggest, as an alternative, a novel way of performing LDA-SIC calculations, which is totally free from self interaction.

## Exchange and correlation in self interaction corrections

The self interaction corrected XC energy, within the approach by Perdew and Zunger, applied to LDA, is normally written as

$$E_{\mathrm{XC}}^{\mathrm{SIC}}[n]=E_{\mathrm{XC}}^{\mathrm{LSD}}[n]\underset{i\sigma }{\overset{\mathrm{occ}}{}}\left[\frac{1}{2}d^3rd^3r^{}\frac{\rho _i^\sigma (𝐫)\rho _i^\sigma (𝐫^{})}{𝐫𝐫^{}}+\frac{1}{2}d^3rϵ_{\mathrm{XC}}(\rho _i^\sigma (𝐫),0)\rho _i^\sigma (𝐫)\right].$$
(2)
The basic idea behind this equation is to subtract the interaction energy of an electron with itself from the total energy. The second term on the right hand side of Eqn.(2) corrects for the self-Hartree interaction but we choose to include this correction to the XC functional since formally one may introduce all terms there which are not found in the kinetic energy functional and the Hartree functional. The last term on the right hand side of Eqn.(2) is called the self XC. The 0 in $`ϵ_{\mathrm{XC}}(\rho _i^\sigma (𝐫),0)`$ indicates that the correction is made separately for each spin-channel. Minimization of the SI corrected energy functional with respect to the density leads to Kohn-Sham like equations with a SI corrected effective potential, $`V_{\mathrm{eff}}`$. The corrections coming from SI are normally referred to as the SIC part of the potential. The correction for the self Hartree interaction is exact, but as will be argued below the correction for self XC is imperfect. To illuminate the problem, we consider an LDA functional (Eqn.(1)) and take a closer look at the XC functional. It is written as $`E_{\mathrm{XC}}(\rho (𝐫),0)=d^3rϵ_{\mathrm{XC}}(\rho (𝐫),0)\rho (𝐫)`$, where one term in the XC energy density, $`ϵ_{\mathrm{XC}}`$, is proportional to $`\rho ^{1/3}`$. Let us now write the total density as $`\rho _{\mathrm{tot}}=\overline{\rho }+\rho _i^\sigma `$, where $`\rho _i`$ is the density given by the orbital which we are interested to correct for self interaction, $`\overline{\rho }`$ is the summed density of all other orbitals and $`\rho _{\mathrm{tot}}`$ the total electron density. It is now easy to see that the correction of PZ, as defined in Eqn.(2), does not remove all self XC, since with the above expressions for $`\rho _{\mathrm{tot}},\overline{\rho }`$ and $`\rho _i^\sigma `$, we obtain from Eqn.(2) an expression for the exchange part of the energy functional

$$d^3r(\overline{\rho }+\rho _i^\sigma )^{4/3}d^3r(\rho _i^\sigma )^{4/3}.$$
(3)
Eqn.(3) is not equal to $`d^3r(\overline{\rho })^{4/3}`$ and hence does not represent a contribution to the energy functional that is free from self interaction.

In a similar way the correction of PZ does not remove the self XC interaction in the effective potential of the Kohn-Sham equations. The XC part of the effective potential, $`V_{XC}`$, of the Kohn-Sham equations is in LDA written as

$`V_{\mathrm{XC}}^{\mathrm{LSDA}}(\rho _i,𝐫_s,x)=\left[\mu _x^P(𝐫_s)+v_c(𝐫_s)\right]2\rho _i^{1/3}+\mu _c^P(𝐫_s)v_c(𝐫_s),`$ (4)
where it may be seen that the exchange potential ($`\mu _x^P(𝐫_s)2\rho _i^{1/3}`$ contains a term proportional to $`\rho (𝐫)^{1/3}`$ (for a description of terms in Eqn(4),see Ref. ). The PZ corrected effective potential for electron state $`i`$ involves, analogous to Eqn.(2), the following terms,

$$V_{\mathrm{eff},i}^{\mathrm{PZ}}=\frac{\overline{\rho }(𝐫^{})}{𝐫𝐫^{}}d^3r^{}+\mathrm{A}(𝐫)([\overline{\rho }(𝐫)+\rho _i^\sigma (𝐫)]^{1/3}\rho _i^\sigma (𝐫)^{1/3}),$$
(5)
and we observe that due to the non-linear dependence on the electron density, $`V_{\mathrm{eff},i}^{\mathrm{PZ}}`$ is dependent on $`\rho _i`$ and hence SI effects are present in the effective potentials well. This was also pointed out by other authors.

## A novel method of self interaction

We have now pointed out why the SIC form of PZ does not remove all the SI and below we suggest a method which does this. In general, the approach is to identify all contributions in the LDA energy expression (Eqn.(1)) and the effective potential which contain SI and to remove these interactions. Writing $`\rho _{\mathrm{tot}}=\overline{\rho }+\rho _i^\sigma `$, one may obtain an energy expression which is free from self-interaction, by applying the correction to the densities instead of the energies. The energies for $`\rho _i^\sigma `$ are calculated from $`\overline{\rho }`$ instead of the total density $`\rho _{\mathrm{tot}}`$. LDA contains a term $`E[\overline{\rho }+\rho _i^\sigma ](\overline{\rho }+\rho _i^\sigma )`$ whereas the wanted expression is $`E[\overline{\rho }+\rho _i^\sigma ]\overline{\rho }+E(\overline{\rho })\rho _i^\sigma `$ ($`E`$ denoting the total energy functional) if we only apply SIC to orbital $`\rho _i^\sigma `$. The expression for the total energy, when formulated as $`E^{\mathrm{tot}}=E^{\mathrm{LDA}}+E^{\mathrm{corr}.}`$, is

$`E_{\mathrm{SIC}}^{\mathrm{tot}}=E_{\mathrm{LDA}}^{\mathrm{tot}}(ϵ^{\mathrm{SIC}}){\displaystyle \frac{1}{2}}{\displaystyle \underset{i\sigma }{}}{\displaystyle \frac{\rho _i^\sigma (𝐫)\rho _i^\sigma (𝐫^{})}{𝐫𝐫^{}}d^3rd^3r^{}}`$ (6)
$`{\displaystyle \underset{i\sigma }{}}{\displaystyle [ϵ_{\mathrm{XC}}(\rho )ϵ_{\mathrm{XC}}(\overline{\rho })]\rho _i^\sigma d^3r}+{\displaystyle \underset{i\sigma }{}}{\displaystyle [V_{\mathrm{eff}}(\rho )V_{\mathrm{eff}}(\overline{\rho })]\rho _i^\sigma d^3r},`$ (7)
where $`E_{\mathrm{LDA}}^{\mathrm{tot}}(ϵ^{\mathrm{SIC}})`$ is the normal LDA total energy expression (Eqn.(1)) evaluated with a KS-eigenvalue obtained from Eqn.(9) below, the effective potential. To arrive at the expression in Eqn.(6) we have explicitly subtracted the parts in the LDA expression which contain SI and thereafter added a SI free part.

The effective potential to be used in the KS-equation is constructed in the usual way

$$V_{\mathrm{eff}}(\rho ,\rho _i^\sigma ,𝐫)=\frac{\delta [E_\mathrm{H}+E_{\mathrm{XC}}]}{\delta \rho _i^\sigma (𝐫)}+V_{\mathrm{eN}},$$
(8)
where $`V_{\mathrm{eN}}`$ is the electron-nucleus potential 2Z/r. Explicitly taking the functional derivative, we get the potential

$`V_{\mathrm{eff}}(\rho _{tot},\rho _i,𝐫)={\displaystyle \frac{\rho _{tot}(𝐫^{})\rho _i^\sigma (𝐫^{})}{𝐫𝐫^{}}d^3𝐫^{}}+\mathrm{A}(𝐫)[\rho _{\mathrm{tot}}(𝐫^{})\rho _i^\sigma (𝐫^{})]^{1/3}+V_{\mathrm{eN}}=`$ (9)
$`{\displaystyle \frac{\overline{\rho }(𝐫^{})}{𝐫𝐫^{}}d^3𝐫^{}}\mathrm{A}(𝐫)\overline{\rho }(𝐫)^{1/3}+V_{\mathrm{eN}}.`$ (10)
This effective potential is by construction totally free from SI if one uses the electron density $`\overline{\rho }`$ in the functional A$`(𝐫)`$. The Hartree part in this equation coincides exactly with the Hartree part in Eqn.(2). In short, the potential for the orbital $`\rho _i^\sigma `$ can be written as $`V(\rho _{\mathrm{tot}}\rho _i^\sigma ,𝐫)=V(\overline{\rho },𝐫)`$, instead of the normal orbital independent potential $`V(\rho _{\mathrm{tot}},𝐫)`$, or the PZ-SIC: $`V(\rho _{\mathrm{tot}},𝐫)V(\rho _i^\sigma ,𝐫)`$. This procedure is applicable to all approximations for the XC energy but here, for simplicity, we choose to discuss the LDA. In the asymptotic limit $`𝐫\mathrm{}`$ the effective potential has the correct $`\frac{1}{r}`$ dependence, since we subtract one electron from the screening cloud of the neutral atom. One can arrive at Eqn.(9) by starting from the Hartree-Fock expression for the effective potential, which is free from self interaction, and then replace the Fock part by A$`(𝐫)\overline{\rho }(𝐫)^{1/3}`$. This replacement then hopefully both simplifies the effective potential by making it local as well as includes correlation effects.

## Numerical tests on the method

We have implemented this new method, and the method by PZ, to be able to make a comparison between the two, in a full potential linear muffin-tin orbital (FP-LMTO) method. Since Pr is a material with 2 localized f-electrons which has been studied extensively by LDA, PZ-SIC and the orbital polarization method we choose to test our theory for this element. For simplicity Pr was put in an fcc structure (observed structure is the dhcp structure), and we used 2197 k-points in the whole Brillouin zone. The SI corrected functionals were implemented only for the 4f core electrons. In the implementation we have ignored the non-orthogonality of the orbitals introduced by SIC since in atomic like problems it has been shown to be negligible. The results for praseodymium are shown in the TableI where the total energy contributions are listed, and in Figures 1 and 2 where the electron density and different contributions to the effective potential are displayed.

By inspecting the 4$`f`$-density (shown in Fig.1) one observes that the orbitals become more contracted when SIC is used. Notable is also that the radial extent of the orbital in the new SIC method is largely reduced compared to the PZ-SIC method. We also note that compared to an atomic Hartree-Fock density the $`f`$-cloud is actually even more contracted. From Fig.2b one sees that for the current SIC method the main contribution to the potential entering the KS equation comes from the Hartree part, the XC part provides a very small contribution. This can be compared with the PZ-SIC potential shown in Fig.2a, where the same Hartree part is used, but there is a big correction to the XC potential. The reason why the correction to the XC potential in this implementation is so small is that the potential has the form $`V_{\mathrm{XC}}=V_{\mathrm{XC}}(\rho _{\mathrm{tot}}\rho _i^\sigma )`$ instead of the LDA $`V_{\mathrm{XC}}=V_{\mathrm{XC}}(\rho _{\mathrm{tot}})`$ and for large values of $`\rho _{\mathrm{tot}}`$ the difference between them, $`\rho _{\mathrm{tot}}\rho _i^\sigma =\overline{\rho }\rho _{\mathrm{tot}}`$, is very small. As a side remark we note that the electronic structure and density of states of the valence states is almost unchanged (no SIC is applied to the valent electrons), showing that the contraction of the 4$`f`$cloud due to SIC only marginally influences the energy-band dispersion of the valence states.

From the TableI it is seen that the KS eigenvalue is significantly lower in the new SIC method compared to LDA and PZ-SIC. The fact that the $`f`$-orbital localizes more gives rise to a larger screening for the valence electrons, which as a result become more extended and for this reason a larger equilibrium volume is found. However, as is seen from the table the difference in equilibrium volumes between the different approximations is only a few percent and they agree rather well with experiment. The small deviation from experiment can be explained by a contribution to the cohesive energy from delocalized $`f`$-electrons as described previously. The cohesive energy was previously calculated from LDA to be 5.09 eV/atom, a value which is 0.55 eV/atom higher than experiment. Our calculation lowers the cohesive energy with 0.29 eV/atom compared to the LDA calculation, and is hence in better agreement with experiment. The PZ calculation lowers the cohesive energy by only with 0.17 eV/atom compared to LDA.

In the Table we also note that the self-Hartree part is larger in our method of SIC than in the PZ-method, this is due to the more contracted 4$`f`$-orbital. In the PZ-SIC the total energy changes only by 0.25 Ry since the Hartree correction is balanced by a change in the XC energy. In our implementation of SIC the Hartree energy is not balanced to the same extent by the XC energy, thus giving a change in total energy by 1.90 Ry. We also compared the 4$`f`$-eigenvalue with an atomic, relativistic Hartree-Fock calculation on Pr3+ with a 4f2 configuration which is free from SI (but without correlation effects). This calculation gives a value of $`ϵ_{4f}`$=-2.870 Ry, which may be compared to an atomic-like calculation (a band-calculation with lattice expanded) which gives -2.3062 Ry for the new SIC method, -1.175 Ry for the PZ-SIC method and -0.444 Ry for a normal LDA calculation. Hence the current SIC method is closer to a Hartree Fock description than LDA or PZ-SIC. This finding is consistent with the discussion above that Eqn.(9) may be viewed as an approximation of the Hartree-Fock equation.

The present SIC method lowers the total energy compared to the method of Perdew and Zunger. For this reason it may not be as suitable for describing localization-delocalization transitions in the rare-earth metals, since the PZ-SIC method reproduces the Mott transition in Ce very well. However, since the PZ-SIC method does not remove all self interaction it would seem that the good description of this method in describing the localization-delocalization transitions in the rare-earths is somewhat fortuitous.

We also remark that as long as a Hartree term (which depends linearly on the electron density) is introduced in the energy functional, any local XC functional $`ϵ_{\mathrm{XC}}(\rho _i(𝐫),0)`$ which is a functional of the total density must also have a linear dependence on the electron density, if SI is to be absent in the functional.

By construction the LDA XC potential cancels some of the unphysical self interaction, e.g. in hydrogen 95$`\%`$ of the self interaction is cancelled. Application of the presented SIC method to the hydrogen-atom would yield an exact result (as would the PZ-SIC). Since the application of PZ-SIC normally improves the LDA result, we believe that the currently described method of SIC will be successful when applied to systems like the 1-D Hubbard model and the bandgap problem in NiO, since it properly corrects for all self interaction. It may also be a good starting point for finding an exact XC functional of localized or nearly localized electron systems.

## Conclusion

In conclusion, we have presented a novel way to perform SIC-corrections in DFT calculations, and compared it to an older, not fully SI free, method. In its form the suggested correction is similar to the corrections of the non-spherical 4f density proposed by Brooks et al. in calculations of crystal field excitations. The method presented is closer to the suggestion of Nesbet and differs from the method of e. g. Manoli et al whose approximation is more similar to that of Ref.. We have applied the new SIC method to praseodymium metal and showed that it results in stronger localization of the f-shell, thus correcting more for self interaction. Despite their large differences in describing the 4f electrons, the calculated equilibrium volumes of LDA, PZ-SIC and the presently proposed SIC method are all quite similar and agree with the experimental values. This is a result of that the valence states (that are not calculated using any SIC potential) are rather independent of the shape of the 4$`f`$-density and it is namely the valence states that are responsible for the chemical bonding. On the other hand, the cohesive energy is better described by the new SIC method, both compared to LDA and to PZ-SIC.

## ACKNOWLEDGMENT

This project has been financed by the Swedish Natural Science Research Council (NFR). We are also thankful to Sverker Edvardson Mitthögskolan, Sundsvall for providing results from the Hartree-Fock calculations.