text:
# Baryons in the Warm-Hot Intergalactic Medium

## 1. Introduction

Observations indicate that most of the baryonic matter in the universe does not reside in galaxies. At high redshifts ($`z2`$), the overwhelming majority of baryons are in a diffuse, photoionized intergalactic medium (IGM), observable as H I absorption lines in the spectra of distant quasars (Cen et al. (1994); Zhang, Anninos, & Norman (1995); Miralda-Escudé et al. (1996); Hernquist et al. (1996)). The baryonic density inferred from these observations (Rauch et al. (1997)) is in good agreement with nucleosynthesis arguments based on observed deuterium abundances (Tytler, Fan, & Burles (1996)). But by redshift zero the total baryonic component, inferred from H I absorption, gas and stars in galaxies, and other observations, has declined to a level small compared to that seen at high redshift and expectations from nucleosynthesis (Fukugita, Hogan & Peebles (1998); Hogan (1999)). Thus the question arises, where are the baryons at the present epoch?

By the current epoch, hierarchical structure formation has produced deep potential wells into which the baryons accrete, thereby moving a significant portion of the baryons from the IGM into stars, galaxies, groups, and clusters. These complex evolutionary processes may now be modeled directly using cosmological hydrodynamic simulations, enabling an investigation into the location and phase of baryonic constituents in the present-day universe, and suggesting possible avenues for their direct detection.

These hydrodynamic simulations of structure formation indicate that baryons in the universe reside in four broad phases, defined by their overdensity $`\delta \rho /\overline{\rho }1`$ (where $`\overline{\rho }`$ is the mean baryonic density) and temperature $`T`$:

1. Diffuse: $`\delta <1000`$, $`T<10^5`$ K. Photoionized intergalactic gas that gives rise to Lyman alpha absorption.
2. Condensed: $`\delta >1000`$, $`T<10^5`$ K. Stars and cool galactic gas.
3. Hot: $`T>10^7`$ K. Gas in galaxy clusters and large groups.
4. Warm-Hot: $`10^5<T<10^7`$ K. The “Warm-Hot Intergalactic Medium” (WHIM), discussed here.

.

Cen & Ostriker (1999; hereafter CO (99)) and Davé et al. (1999; hereafter DHKW ) predicted that a sizeable fraction of all baryons at the present epoch reside in this last warm-hot phase (see Figure 2 of CO (99) and Figure 12 of DHKW ). Such a reservoir has significant implications for producing an accurate census of baryons for comparison with nucleosynthesis arguments (e.g. Fukugita, Hogan & Peebles (1998)) because gas at these temperatures and densities is difficult to detect in either absorption or emission, as we will discuss in §6.

In this paper we study the nature and evolution of warm-hot gas in a $`\mathrm{\Lambda }`$CDM universe, using a suite of cosmological hydrodynamic simulations having a wide range of numerical and physical parameters. The purpose of this paper is to ask how robust simulation predictions are to these parameters, to examine the physical state of warm-hot gas in the universe, and to investigate constraints on warm-hot gas from soft X-ray background observations.

In §2 we briefly describe the simulations used. The evolution of diffuse gas at high redshift into condensed, hot, and WHIM gas at the present epoch is quantified in §3. Our primary result, presented in §4, is that the WHIM accounts for a significant fraction ($`3040\%`$) of baryonic mass at $`z=0`$, regardless of variations in spatial resolution, input physics, or hydrodynamic algorithm. This is because the evolution of the WHIM is driven primarily by shock-heating of gas falling into gravitationally-generated potential wells, a process which is well-understood and modeled, and only secondarily by processes such as supernova feedback, radiative cooling, and photoionization. Furthermore, most of the WHIM gas is at relatively low overdensities, so shock heating of intergalactic gas occurs during flows onto non-equilibrium large-scale structures such as filaments. The majority of warm-hot gas is found outside of virialized objects such as galactic halos and galaxy groups. The low overdensities explain why radiative cooling and supernova heating do not drive its evolution, and why the presence of this component probably does not violate constraints from the X-ray background, as we show in §5. Finally, we briefly discuss strategies for direct detection of this gas, noting that the easiest place to detect emission from WHIM gas is relatively close to galaxies, where it is dense, even though most of the warm-hot baryons are not in these regions. In summary, the WHIM is a robust and generic prediction of currently popular hierarchical structure formation models, and it contains roughly one-third of all baryons in the universe today.

## 2. Simulations

We use six cosmological hydrodynamic simulations of randomly-selected volumes in $`\mathrm{\Lambda }`$-dominated Cold Dark Matter universes, employing three different numerical techniques, with a range of physical and numerical parameters. These parameters are summarized in Table 1. Simulations run with Parallel TreeSPH (Davé, Dubinski & Hernquist (1997)) are labeled D1 and D2, simulations run with TVD-PM (Ryu et al. (1993)) are labeled C1 and C2, and simulations using Adaptive Mesh Refinement (AMR; Bryan (1999)) are labeled B1 and B2. The cosmology chosen is close to the “concordance model” which is in agreement with a wide variety of observations (Bahcall et al. (1999)).

For our purposes, the significant inputs are the variations in spatial resolution ($`1200h^1\mathrm{kpc}`$), baryonic mass resolution ($`2\times 10^710^9M_{}`$), boxsize ($`11.11100h^1\mathrm{Mpc}`$), input physics, and hydrodynamic algorithms. D1 and D2 are high spatial resolution Lagrangian (particle-based) simulations, C1 and C2 are lower spatial resolution Eulerian simulations having high mass resolution and employing the Total Variation Diminishing scheme, and B1 and B2 are high resolution adaptive mesh simulations based on the Piecewise Parabolic Method (note the exceptional resolution of $`1h^1\mathrm{kpc}`$ achieved by AMR in the high density regions). All simulations except B1 include radiative cooling from H and He, photoionization heating, and star formation; B1 includes none of these. C1, C2, and B2 additionally include metal-line cooling, with the metallicity determined self-consistently from supernova output. All simulations use a $`\mathrm{\Lambda }`$-dominated cold dark matter universe, having similar power at cluster and galaxy scales. B2 has only been evolved to $`z=0.75`$, but this will be sufficient to indicate the relevant trends.

## 3. Evolution of Intergalactic Gas

Figure 1 shows the evolution of the baryonic mass fraction in the four phases described above. The four panels show results from simulations D1, D2, C1 and C2. Despite differences in simulation volume, resolution, and numerical method, the evolution of various phases is qualitatively similar. At high redshift ($`z2`$), the dominant fraction of baryons resides in diffuse gas (dashed lines) giving rise to Lyman alpha forest absorbers, as has been explored in detail elsewhere (see Rauch (1998) for a review). As structure forms, diffuse gas is shock-heated, producing warm-hot gas (solid lines). Gas that is driven to higher densities due to gravitational instability is able to cool into the condensed phase (dotted lines) and form stars. At lower redshifts ($`z2`$), large potential wells are produced that shock-heat gas to $`T>10^7`$ K, giving rise to hot cluster gas (dot-dashed curves). Nevertheless, at the present epoch, the total fraction of baryons in clusters is small.

Figure 1: Evolution of mass fractions in four baryonic phases, in four simulations.

Figure 1 shows that the fraction of diffuse baryons at the present epoch is between 20% and 40%, the fraction of warm-hot baryons is $`3040\%`$. Gravitationally bound baryons, i.e. those in stars, galactic gas, clusters, and intragroup media, make up the rest. Thus baryons in the present-day universe are divided roughly equally among diffuse, warm-hot, and bound components.

While all the simulations are qualitatively similar, the exact distribution of gas in these phases is sensitive to details of the simulations. For instance, the condensed phase fraction is sensitive to how gas cools and forms stars in these simulations, which in turn is significantly affected by resolution (since cooling and star formation $`\rho ^2`$). Also, the effect of supernova feedback heating depends strongly on resolution, as we will discuss later. The growth of structure is affected by the amount of large-scale power present. In particular, the hot gas fraction in clusters is sensitive to cosmic variance, since it is dominated by the largest virialized objects in the volume. Other phases, including the warm-hot phase, are less sensitive to volume effects, as we show in §4.1.

All things considered, it is not surprising that there are differences up to a factor of two in the fractions in various phases at $`z=0`$. While the differences may be significant, an investigation of their exact causes is beyond the scope of this work, though these differences offer clues into the physical processes driving WHIM evolution, as we will explore in §4.1. Rather, we focus on the remarkable qualitative consistency in the evolution of various gas phases, given the variety of simulation methodologies utilized.

For the condensed phase, a census of baryons in stars and cold gas estimates its mass fraction to be around 20% (Fukugita, Hogan & Peebles (1998)). The simulations shown generally produce somewhat higher values for the condensed fraction, but there are uncertainties in the contributions from low-mass stars and supernova remnants. Of these simulations, C2 contains the largest fraction of condensed gas ($`40\%`$), while the rest are all lower, down to $`25\%`$ for simulation C1. B1, of course, has virtually no condensed gas since it does not include cooling. In what follows, we will fix the condensed fraction to be 20% and redistribute the excess condensed gas equally among all the other components, assigning a somewhat arbitrary fraction of 1/3 of it to the WHIM (the only component we examine from here on), in order to facilitate a less resolution-dependent comparison of the intergalactic components in these simulations.

## 4. The Warm-Hot Intergalactic Medium

### 4.1. Evolution of $`\mathrm{\Omega }_{\mathrm{WHIM}}`$

Figure 2: Evolution of $`\mathrm{\Omega }_{\mathrm{WHIM}}`$ in our six simulations. For this comparison, the baryonic fraction in cold galactic gas and stars has been fixed at 20% in all simulations.

Figure 2 shows the evolution of the mass fraction of baryons in warm-hot gas, $`\mathrm{\Omega }_{\mathrm{WHIM}}/\mathrm{\Omega }_b`$, for our six simulations, with the condensed component having been fixed at 20%. C1 and C2 have the highest WHIM fractions, at 42% and 35% of baryonic mass, respectively, at $`z=0`$. D1 and D2 have lower WHIM fractions at the present epoch, around $`30\%`$. B2 has only been evolved to $`z=0.75`$, but its evolution closely mirrors that of the Lagrangian runs, so we expect its WHIM fraction will also be $`30\%`$ by $`z=0`$. These values are listed in Table 1. By the present epoch, all simulations have WHIM fractions within 50% of each other. Simulation B1 does not include radiative cooling, so it is an unphysical model that will only be used for comparison with B2.

We expect that the evolution of warm-hot gas is governed by shock heating of intergalactic gas onto large-scale structure, supernova feedback, and radiative cooling. Numerical considerations such as resolution, volume and hydrodynamic algorithm may also play a role. In this section we use case-by-case comparisons among our six simulations to examine how each of these processes affects WHIM evolution in our models.

All of our simulations explicitly include the growth of structures and the accretion of gas onto those structures. The qualitative consistency of WHIM evolution in all our simulations, including simulation B1 without cooling or star formation, suggests that gravitational shock heating of gas falling on large-scale structures is the dominant heating mechanism for WHIM gas. This is an important point, as it suggests that all other differences between these simulations are of secondary importance.

The most obvious differences between models are that D1, D2 and B2 all predict similar WHIM evolution, while C1 and C2 have more WHIM gas. The primary distinction between these sets of simulations is spatial resolution. This manifests itself in various ways, as it affects gas cooling, the rate at which stars form, and the the injection of supernova heat energy into intergalactic gas. To elaborate on this last point, high-resolution runs (D1, D2, B2) deposit feedback energy (thermal only) locally in very high density regions where stars are forming, and thus it quickly radiates away, with virtually none of it being distributed into the intergalactic medium. Conversely, lower-resolution runs (C1 and C2) deposit the same supernova energy over hundreds of kiloparsecs, resulting in a significant fraction of feedback energy escaping into diffuse regions where it cannot radiatively cool away. None of the simulations here are capable of resolving supernova-driven galactic winds through a multi-phase interstellar medium, which are likely to be responsible for distributing energy and metals into the diffuse IGM (Mac Low (2000); Efstathiou (2000)), thus we are relying on heuristic modeling of these processes. As we discuss in §4.2, it is in the distribution of feedback energy where resolution plays its most crucial role in predicting the evolution of WHIM gas.

First we examine the effect of radiative cooling. Simulation B1 has no cooling, whereas simulation B2 is a similar run with cooling (and is our highest spatial resolution run). Cooling has a greater effect at earlier times because intergalactic gas is denser then and thus can cool significantly. After $`z3`$, the rate of growth of the WHIM fraction is similar in B1 and B2, suggesting that WHIM gas is no longer affected by cooling. Even by $`z=0.75`$, the difference between B1 and B2 is not large, indicating that cooling plays a minor role in the overall evolution of WHIM gas.

Eulerian codes resolve shock fronts better than SPH in diffuse regions (Kang et al. (1994)), because they have a higher density of resolution elements there, and, in the case of the Eulerian codes discussed here, because they incorporate explicit shock capturing algorithms to resolve fronts over two cells. If a significant component of shock heating arises from small-scale shocks unresolved by SPH, Eulerian codes may produce higher temperatures. One might suspect that the difference between PTreeSPH runs and TVD-PM runs could be partially due to this effect. However, B2 is also an Eulerian simulation, yet it has nearly the same WHIM fraction as the Lagrangian runs. So this cannot be a significant effect.

Simulation volume could also play a role in the amount of WHIM gas, since larger volumes contain larger perturbations that can result in stronger gravitational shock heating. This may explain the difference between C1 and C2, whose box lengths differ by a factor of two. However, D1 has a box length five times that of D2, yet their WHIM fractions are nearly identical at all times. Furthermore, B2 has twice the box length of D1, yet its WHIM fraction is in good agreement with D1 and D2. As we will show in the next section, typical WHIM gas is at very moderate overdensities, and cosmic variance in that regime is typically small, in contrast to the hot IGM fraction which is dominated by rare, massive objects (i.e. clusters) whose numbers are more sensitive to simulation volume.

It is possible that opposing effects are causing D1, D2 and B2 to be similar. For instance, D1 resolves shock fronts better than D2, perhaps resulting in more shock heating and making up for its lack of large-scale power. It is also possible that the adaptive refinement of B2 results in more cooling in shock fronts, compensating for its larger volume. However, unless all of these effects are of secondary importance, cancellations at the level we find here would require a remarkable coincidence.

The parameters of the underlying cosmological model are expected to have a non-negligible effect on the WHIM component. The rate of structure evolution in large part determines how much shock-heated gas is present at any epoch. While we have not sampled a range of cosmologies here, DHKW  examined four cosmologies, namely $`\mathrm{\Lambda }`$CDM, Tilted CDM, Cold+Hot DM, and Open CDM models, using numerical parameters similar to those of simulation D2. The differences between the WHIM fractions (called “shocked” gas in DHKW ) in those models at the current epoch is comparable to the differences seen here due to other factors, as indicated by their Figure 12. Similarly, CO (99) examined $`\mathrm{\Lambda }`$CDM, Cold+Hot DM, and Open CDM models, and also found broad consistency. Thus while our quoted fractions may be specifically for a $`\mathrm{\Lambda }`$CDM cosmology, our qualitative conclusions are unlikely to be highly sensitive to cosmology.

In summary, the evolution of $`\mathrm{\Omega }_{\mathrm{WHIM}}`$ is qualitatively consistent among all simulations examined, and results in $`3040\%`$ of baryons residing in the WHIM today. Radiative cooling, simulation volume, and algorithmic details do not significantly affect WHIM evolution in these runs.

Figure 3: WHIM gas in simulation C2. Contours are color coded by overdensity; green represents overdensity $`\delta 10`$, while red shows $`\delta 10^4`$.

### 4.2. The Physics of WHIM Gas

We now explore the physics that drives the formation and evolution of WHIM gas. A qualitative physical picture of WHIM gas may be obtained by examining Figure 3, which shows the location of gas in simulation C2 having $`10^5<T<10^7`$K, with contours color-coded by density. WHIM gas is seen to primarily trace out filamentary large-scale structures. Like the intergalactic medium, WHIM gas does cluster around dense regions that are sites of galaxy formation. However, we will show below that the majority of WHIM gas is contained in the filaments.

Figure 4: Top panel shows the mass fraction of WHIM gas as a function of density at $`z=0`$ in our simulations. Bottom panel shows the same quantity at $`z=0,1,2,3`$ in simulation D1.

The top panel of Figure 4 quantifies the spatial distribution of WHIM gas in the universe. It shows a histogram of WHIM gas mass as a function of density, for simulations D1, D2, C1, and C2. All four simulations consistently show that the dominant fraction of WHIM gas is at relatively low densities, with a peak around an overdensity of $`1030`$. Simulations B1 and B2 are not shown because B1 does not include cooling and B2 has only been evolved to $`z=0.75`$. Still, B1 at $`z=0`$ and B2 at $`z=0.75`$ show peak overdensities of 18 and 21, respectively, so they are consistent with the other simulations. $`7080\%`$ of WHIM baryons lie in the overdensity range $`5<\delta <200`$, typical of filaments.

The typical overdensities of WHIM gas are much smaller than that of matter contained in bound, virialized objects. A maximal estimate of the bound fraction of WHIM gas may be obtained as the fraction of WHIM gas with $`\delta 60`$, which is approximately the overdensity at the virial radius of an isothermal sphere in a $`\mathrm{\Lambda }`$CDM cosmology. The mass fraction of WHIM gas with $`\delta 60`$ is $`30\%`$ in our simulations. Clearly some of this gas will not be bound, but rather infalling material. As an independent check, we use the group finding algorithm SKID111http://www-hpcc.astro.washington.edu/TSEGA/tools/skid.html (Spline Kernel Interpolative DENMAX) to identify bound warm-hot particles in simulation D1, and find a bound fraction of WHIM gas between $`10\%`$ and $`25\%`$, depending on the linking length used ($`50500h^1\mathrm{kpc}`$). The lower end of this range probably corresponds to gas contained in galactic halos, although the extent of a galactic halo becomes ill-defined when it resides within a group or cluster. In summary, WHIM gas is mostly an intergalactic component, with a majority of it residing outside of virialized structures such as galaxies or groups. Note that coronal gas in galaxies also lies in the warm-hot temperature range, but it is a very small fraction of the galactic baryonic mass, as most galactic baryons are tied up in stars and cold gas.

How can intergalactic gas be heated to $`T>10^5`$ K without being bound in a massive virialized halo? One way would be if energy was added from non-thermal processes such as supernova feedback. However, this is not the driving process for WHIM gas in these simulations. Supernova feedback energy is added in all simulations shown in Figure 4, but in simulations D1, D2 and B2 the effect is negligible, since, as discussed before, feedback is added purely thermally into very dense regions, where it radiates away almost immediately and adds no heat to diffuse regions. But even in these high-resolution simulations, where supernovae add negligible heat to WHIM gas, the typical overdensities are smaller than those typical of virialized halos.

Supernova heating is, however, likely to be responsible for the somewhat lower overdensities in simulations C1 and C2 as compared to the higher resolution simulations, so its effect is non-negligible. The differences cannot be due to other spatial resolution effects such as cooling, as C1 and C2 are themselves quite similar, as are D1 and D2. Instead, the lower spatial resolution in C1 and C2 results in significant supernova energy being deposited in intergalactic gas, as described in the previous section. This additional feedback heating raises the pressure of intergalactic gas and lowers the typical density. Furthermore, because supernova energy is not radiated away immediately, there is a small component of warm-hot gas in and around galaxies seen as the high-density bump for C1 and C2 in Figure 4. Note that the densities in D1 and D2 are computed slightly differently than in C1 and C2; in the former, Lagrangian runs, the density field is smoothed by the (variable) smoothing length of the particle, while in the Eulerian runs the smoothing is fixed at the cell size. However, since in both cases the smoothing is done on scales smaller than the density variations, this is not expected to produce any systematic differences in the resulting densities.

The effect of supernova heating can be roughly estimated by the following simple argument. Consider two extreme cases, one where all supernova energy is deposited within dense regions where it immediately radiates away, and another where supernova energy is distributed uniformly over all baryons in the universe. The former case is a reasonable approximation for the Lagrangian runs D1 and D2 as well as the high-resolution AMR run B2, while the latter is closer to the Eulerian runs C1 and C2 (though clearly a much more extreme case). In the former case, supernovae add no heat to WHIM gas. In the latter, the specific heat added per baryon yields a temperature increase of

$$\delta T_{\mathrm{SN}}=\frac{\mathrm{\Omega }_{}\rho _{\mathrm{crit}}ϵ_{\mathrm{SN}}}{k_B\overline{n}_H},$$
(1)
where $`\mathrm{\Omega }_{}`$ is the cosmic mass fraction in stars, $`\rho _{\mathrm{crit}}`$ is the critical density, $`ϵ_{\mathrm{SN}}`$ is the specific energy output per unit mass of stars formed, and $`\overline{n}_H`$ is the mean number density of H atoms. We take $`\mathrm{\Omega }_{}=0.1\mathrm{\Omega }_b=0.002h^2`$, and $`ϵ_{\mathrm{SN}}=2.5\times 10^{48}`$ erg $`\mathrm{g}^1M_{}`$ from a Salpeter IMF with each supernova from a star with $`M>8M_{}`$ outputting $`10^{51}`$ ergs. Then, $`\delta T_{\mathrm{SN}}=0.2`$ keV $`2\times 10^6`$ K. Assuming isentropic heat distribution, $`\rho T^{\frac{3}{2}}=`$ constant, so gas with $`T4\times 10^6`$ K (the maximum of the temperature distribution, as we will show in Figure 5) will have its density reduced by roughly a factor of two compared to the no heating case. This is roughly the level of reduction seen in the Eulerian as compared to the Lagrangian runs. Clearly this model is overly simplistic, as it would predict that all intergalactic gas has $`T10^6`$ K, but it does roughly indicate the magnitude of the effect of supernova heating.

Note that there are no other heating processes in these simulations that contribute significantly to the WHIM. For instance, photoionization only heats gas to $`10^4`$ K. Any other heating process one could devise, such as cosmic rays or supernovae occuring in dwarf galaxies in voids, would add more pressure support to the gas, and thus would push the typical density of WHIM gas even lower. Thus our simulations indicate that WHIM gas is heated to $`T>10^5`$ K primarily by shock-heating of gas accreting onto large-scale structure. These structures, typically filamentary, are not virialized or in dynamical equilibrium.

In the previous section we showed empirically that radiative cooling has a minor effect on WHIM evolution. This may be understood physically given that warm-hot gas is typically at such low overdensities that it cannot collapse into virialized objects, and due to the metagalactic photoionizing background most intergalactic gas lies well within the optically thin regime (see DHKW , Figure 10). Thus typical WHIM gas is too diffuse to self-gravitate or self-shield, and has no way to achieve the densities required to make radiative processes important. Improving the numerical resolution of our simulations would not change this result, as it is based on simple physical arguments.

The bottom panel of Figure 4 shows a histogram of WHIM gas mass for the D1 model at $`z=0,1,2,3`$ (solid, dotted, dashed, long dashed lines). There is a slight trend to higher overdensities at earlier times, since at early times the contribution from virialized structures is greater. This is because gas has not had time to accrete and shock on more diffuse structures, and the largest structures at earlier times have lower temperatures that can fall into the warm-hot range. Still, this is a minor effect; basically, the peak overdensity does not evolve significantly with redshift.

Figure 5: The mass fraction of baryons as a function of temperature in simulation D1, at $`z=0`$ (solid), $`z=1`$ (dotted), $`z=2`$ (short dashed) and $`z=3`$ (long dashed). The arrows indicate the predicted peak temperature from gravitational shock heating at various $`z`$, from equation 2.

Figure 5 shows the temperature distribution of intergalactic gas with $`T>10^{4.5}`$ K in simulation D1 at $`z=0,1,2,3`$. The amount of gas at these warm-hot and hot temperatures grows with time, as seen from Figures 1 and 2. The temperature at the peak of the distribution grows in time as well, reflecting the fact that the universe contains hotter structures at later times. The peak temperature at all redshifts falls within the warm-hot range, and by $`z=0`$ it is up to $`4\times 10^6`$ K. Note that this is close to the temperature of the excess diffuse emission seen by Wang & McCray (1993) in ROSAT data. Figure 5 also shows that the mass of the WHIM component is insensitive to our somewhat arbitrary choice of $`10^5`$ K as the defining lower temperature. A choice of $`10^{4.5}`$ K (as in DHKW ) or even $`10^{5.5}`$ K would not drastically affect our conclusions.

An analytic estimate of the evolution of the peak temperature can be obtained by considering the temperature of intergalactic gas shock heated on mildly nonlinear large-scale structure. If the length scale going nonlinear at a given epoch is $`L_{\mathrm{nl}}`$, and the perturbation collapses on a timescale $`t`$, the resulting sound speed behind the shock will be $`L_{\mathrm{nl}}/t`$. Given that the perturbation has taken a Hubble time to collapse, $`tH^1`$, where $`H`$ is the Hubble constant at time $`t`$. The resulting post-shock temperature is then (CO (99))

$$T_{\mathrm{nl}}c_{\mathrm{nl}}^2=K(HL_{\mathrm{nl}})^2,$$
(2)
where $`K`$ is a constant, at any given epoch. The value of $`T_{\mathrm{nl}}`$, with $`K=0.3`$, is shown at each redshift plotted in Figure 5 by the arrows above each curve. This value of $`K`$ produces roughly the correct peak temperature at $`z=0`$, and as can be seen from Figure 5, it is a reasonable fit to the peak temperature at higher redshifts, though it becomes progressively more difficult to identify a peak. This value of $`K`$ also produces the correct evolution of the globally averaged temperature in simulations C1 and C2, as shown in CO (99). Thus the evolution of the temperature distribution of WHIM gas is consistent with the interpretation that it is heated by gravitationally-induced shocks on mildly nonlinear large-scale structure. Since there is a wide range in the properties of the collapsing structures and therefore infall velocities, there is also a wide range in gas temperatures.

Figure 6: Contours in temperature and density for simulation D1 at $`z=0`$, enclosing 10%, 50% and 90% of the baryons in the range shown. Density and temperature are correlated in the WHIM regime. Thick line indicates a scaling of $`\rho /\overline{\rho }_b=T/10^{4.7}`$ in the warm-hot temperature range.

The temperature and density of WHIM gas are correlated. Figure 6 shows a contour plot of mass within the warm-hot range, as a function of density and temperature, for the D1 model at $`z=0`$. The contour levels enclose 10%, 50% and 90% of the mass in the temperature and density ranges shown in the plot. The thick line indicates an “equation of state” $`\rho T`$ that provides a reasonable fit to gas in the range $`10^5<T<10^7`$ K. This relationship is different from that of diffuse gas, which typically has $`\rho T^{1.7}`$, and the temperature-density relation of WHIM gas has much greater scatter. The higher temperature, different slope, and greater scatter all reflect the importance of shock heating as the dominant mechanism controlling the thermal properties of WHIM gas; the “equation of state” for diffuse gas, on the other hand, arises from the competition between photoionization and adiabatic cooling due to Hubble expansion (Hui & Gnedin (1997)). Figure 6 also suggests that detecting WHIM gas in emission will be easier for gas that is at the highest end of the WHIM temperature range, since that gas will be both denser and hotter. However, the dominant portion of WHIM gas is at lower temperatures, which is perhaps most easily detected via absorption lines (Tripp, Savage & Jenkins (2000)).

## 5. Constraints From the Soft X-Ray Background

Gas with temperatures in the range $`10^5<T<10^7`$ K will emit thermally in the soft X-ray band. The extragalactic soft X-ray foreground (SXRB) flux at 0.1-0.4 keV is roughly $`2035\mathrm{keV}\mathrm{cm}^2\mathrm{s}^1\mathrm{sr}^1\mathrm{keV}^1`$ (Warwick & Roberts (1998)), though uncertainties are large because galactic coronal gas provides an increasing foreground to lower energies. At slightly higher energies ($`1`$ keV), the XRB has further been resolved nearly completely ($`8090\%`$) into point sources, mostly AGN (Mushotsky et al. (2000)). Reasonable arguments then allow only a small contribution to the SXRB from diffuse gas, $`4\mathrm{keV}\mathrm{cm}^2\mathrm{s}^1\mathrm{sr}^1\mathrm{keV}^1`$ (Wu, Fabian & Nusser 1999, hereafter WFN ). Such a limit, in principle, places constraints on the amount of gas in the universe at warm-hot temperatures.

These limits were explored in two independent papers using similar methodologies, WFN  and Pen (1999). Both papers argue that the standard picture of hierarchical formation of virialized objects results in a predicted SXRB that exceeds the observed limits. They suggest that significant non-gravitational heating, typically $`1`$ keV per baryon, is required to unbind warm-hot gas from virialized objects in order to satisfy the SXRB constraints. In this section we discuss these constraints in the context of WHIM gas, and find that our simulations paint a very different picture for soft X-ray emission than the simple models assumed in those two papers. A full calculation of the SXRB from these simulations is a complicated undertaking (because of metallicity, bandpass, and numerical issues) that is beyond the scope of this paper. However, the physical properties of WHIM gas in the simulations are quite different from the properties assumed by WFN and Pen (1998), and we will show that scaling our results to theirs suggests that the observed SXRB flux does not rule out the WHIM predicted by these simulations. The essential difference is one of density. Both WFN and Pen (1998) base their calculations on a Press-Schecter  (1974) analysis, which implicitly assumes that gas is in virialized objects with typical overdensity $`200`$. However, most of the WHIM gas in the simulations is in lower-density filamentary structures rather than virialized objects, thus the SXRB emission is lower (cf. Figure 5).

We can quantify this difference in typical overdensity by considering the clumping factor of the emitting gas. If we define the clumping factor for gas component $`g`$ as

$$C_g<\rho _g^2>/<\rho _g>^2,$$
(3)
then the free-free emissivity from that component is

$$ϵ_{\mathrm{SXRB}}<\rho _g^2T_g^{0.5}>C_g\mathrm{\Omega }_g^2T_g^{0.5},$$
(4)
where $`\rho _g`$, $`\mathrm{\Omega }_g`$, $`T_g10^6`$ K, and $`C_g`$ are the density, mass fraction, temperature, and clumping factor of the gas emitting in soft X-rays. The flux $`j_{\mathrm{SXRB}}`$ of soft X-ray background is then the emissitivity multiplied by path length $`\frac{1}{3}cH^1`$. WFN  argue, sensibly, that it is predominantly warm-hot gas ($`10^5<T<10^7`$ K) that is responsible for soft X-ray emission. This means the appropriate clumping factor $`C_g`$ is that of warm-hot gas, $`C_{\mathrm{WHIM}}`$.

There are several ways to calculate $`C_{\mathrm{WHIM}}`$ in our simulations. One can directly calculate it from equation 3, which is the approach we use for our Eulerian simulations. For Lagrangian simulations, because each particle represents a different volume of gas, it becomes more numerically convenient to calculate $`C_{\mathrm{WHIM}}\xi _{\mathrm{WHIM}}(0)`$, where $`\xi _{\mathrm{WHIM}}(r)`$ is the two-point correlation function of WHIM gas at radius $`r`$.

The resulting $`C_{\mathrm{WHIM}}`$ values for our simulations at $`z=0`$ are listed in Table 1. All simulations show clumping factors in the range $`30400`$. The smaller clumping factors in C1 and C2 arise directly because the WHIM gas is typically less dense in these models as compared to D1, D2, B1 and B2 due to greater supernova feedback in the diffuse IGM, as explained before.

Our clumping factors are significantly less than that used by Pen (1999), who adopts $`C900`$ and argues $`C10^4`$. The reason is that the clumping factor of Pen (1999) is actually that of all gas (assumed to trace the dark matter), not the warm-hot baryons. Such a clumping factor is dominated by the contribution from galaxies and other collapsed, virialized structures. This is not the appropriate clumping factor with which to calculate the soft X-ray background, since as we have shown in the previous section, WHIM gas is typically at much lower overdensities than virialized structures. In Pen (1999), the clumping factor is strongly dependent on simulation resolution, because it is dominated by the highest overdensity objects. Amongst our simulations, a wide range of resolutions give barely an order of magnitude difference in clumping factor, since the WHIM gas is predominantly diffuse.

Pen (1999) determines that $`1`$ keV of non-gravitational heating would reduce his clumping factor to $`C60`$, which he argues is necessary to satisfy XRB constraints. However, this was derived assuming that all gas contributes to the SXRB (see his eq. 2), whereas in our simulations only $`30\%`$ of the baryons contribute to the SXRB (probably an even smaller fraction since the narrow range of SXRB energies arises from a narrower range of gas temperatures than $`10^5<T<10^7`$ K). From equation 4, we see that a reduction in $`\mathrm{\Omega }_g`$, the mean density of emitting gas, results in an increase in the maximum allowed clumping factor by $`(\mathrm{\Omega }_b/\mathrm{\Omega }_{\mathrm{WHIM}})^210`$. Thus his analytically-derived constraint translates to $`C_{\mathrm{WHIM}}600`$. Our clumping factors ($`400`$) are consistent with his analysis, despite the fact that WHIM gas in our simulations undergoes almost no non-gravitational heating.

WFN  present a more detailed model from which we can crudely estimate our SXRB by comparison. They do not specifically use clumping factors to predict the SXRB in their model, but they perform a similar calculation to that of Pen (1999), based on a Press-Schecter  analysis with an assumed halo profile. The Press-Schecter  formalism makes two assumptions about how gas is distributed and heated in the universe: (1) All soft X-ray emitting gas is bound in virialized halos (namely, groups of galaxies); (2) Gas can only be heated by accretion onto a virialized halo, with the gas temperature set by the virial temperature of the halo. Unlike Pen (1999), they allow some gas to cool and therefore not emit in X-rays, but this is a small correction. In their model, the gas clumping factor is the clumping factor of virialized objects, since they are assuming that emitting gas solely resides in such objects.

Our simulations suggest that both Press-Schecter  assumptions are strongly violated in the case of WHIM gas. First, only a small fraction of WHIM gas is bound in virialized halos; most is distributed much more diffusely. Second, in our simulations, intergalactic gas can be heated significantly by purely gravitational processes prior to being accreted onto a virialized object.

Figure 7: Two-point correlation functions for WHIM gas from simulation D1 (solid line) and D2 (dashed line), and for galaxies from simulation D1 (dotted line).

The two-point correlation function $`\xi _{\mathrm{WHIM}}`$ illustrates quite clearly that WHIM gas is not bound in virialized halos, and is instead associated with large-scale structure. Figure 7 shows $`\xi _{\mathrm{WHIM}}`$ for simulation D1 (solid line), D2 (dashed), as well as $`\xi _{\mathrm{galaxy}}`$ (dotted line) calculated from D1. $`\xi _{\mathrm{galaxy}}`$ is computed as the correlation function of all star particles in this simulation, where stars are formed slowly out of cold, dense gas and are seen to trace the galaxy population (see Katz, Weinberg & Hernquist (1996)). At scales less than a few hundred kpc, $`\xi _{\mathrm{WHIM}}`$ flattens, indicating that WHIM gas does not cluster on scales smaller than those typical of large-scale structures. Conversely, the galaxy correlation function continues as approximately a power law in to the resolution limit of our simulation. The clumping factor of WHIM gas, $`\xi _{\mathrm{WHIM}}(0)`$, is well-defined since we can meaningfully extrapolate from our resolution limit to $`r=0`$. This is not true for the clumping factor of galaxies, for which we may only set a lower limit (consistent with Pen 1999), $`C_{\mathrm{galaxy}}4\times 10^4`$.

We can crudely estimate the SXRB flux we would predict in comparison to that predicted by WFN  by multiplying their predicted $`j_{\mathrm{SXRB}}`$ by $`C_{\mathrm{WHIM}}/C_{\mathrm{galaxy}}`$. This corrects for the fact that they use virialized objects (which we take as galaxies) instead of the true diffuse WHIM responsible for SXRB emission. Thus we would roughly predict a SXRB flux that is two orders of magnitude lower than that predicted by WFN . Examining their Figure 2 or 3 shows that such a reduction in $`j_{\mathrm{SXRB}}`$ makes the SXRB contribution from warm-hot gas consistent with observational limits.

In summary, the discrepancies in the predicted soft XRB versus previous studies arise from their use of a Press-Schecter  analysis, or an equivalent method that implicitly places emitting gas into virialized objects, in order to study an intergalactic component of baryons that is inherently much more diffuse. These methods do not allow for gravitational shock heating in unbound objects such as filaments, thus they are forced to postulate significant non-gravitational heating to make the emitting gas more diffuse. In our scenario, most WHIM gas has never fallen into virialized objects, though it is still heated (almost) purely via gravitational processes. Scaling previous analyses by the typical overdensities found in our simulations suggests that the properties of the predicted WHIM are consistent with SXRB constraints. However, these scaling arguments are approximate at best, and not enough to guarantee consistency. More accurate calculations of X-ray emission from these simulations are certainly warranted and may yield interesting constraints on the WHIM component; these issues will be addressed in detail in separate papers (Fardal et al., in preparation; Phillips et al., in preparation).

## 6. Summary

We study the warm-hot intergalactic medium (WHIM), defined as all the gas in the universe with temperature $`10^5<T<10^7`$ K, in six cosmological hydrodynamic simulations with widely varying spatial resolutions, volumes, code algorithms, and input physics. In each simulation, the WHIM contains $`3040\%`$ of all baryons in the present-day universe. As a rule of thumb, our simulations predict that the fractions of baryons in the warm-hot phase, the diffuse phase, and gravitationally bound systems are roughly comparable at the present epoch.

The WHIM is comprised primarily of gas at moderate overdensities, with a median overdensity $`1030`$, and is predominantly an intergalactic component. It is shock-heated by accretion onto non-equilibrium filamentary large-scale structures, with possibly a small energy contribution from non-gravitational processes such as supernova feedback. Despite being heated primarily by gravitational processes, WHIM gas in our simulations is consistent with constraints from the soft X-ray background. The clumping factors of WHIM gas in our simulations range from $`30400`$, which are far below that of virialized objects ($`10^4`$) that some studies assume are the sources of soft XRB photons.

Different simulations give somewhat different fractions of WHIM gas at the present epoch, though the fraction is significant in all simulations. We argue that the differences are primarily attributable to resolution effects, specifically manifested in the distribution of supernova feedback energy in different simulations. In high spatial resolution simulations, the added heat remains in dense regions and radiates away almost immediately, whereas in lower spatial resolution simulations a significant fraction of heat ends up in diffuse intergalactic gas where it cannot radiate away. Still, the differences in present-day WHIM fraction are not large between these two somewhat extreme cases for feedback. Radiative cooling plays a minor role in the evolution of warm-hot gas because it is mostly at low densities, and is too diffuse to self-gravitate and self-shield into mini-halos. Other physical effects such as star formation, photoionization, and cosmic variance do not produce significant differences in the amount or properties of warm-hot gas in our simulations.

Since the exact predicted fraction of WHIM gas is as yet sensitive to simulation details, it would be greatly beneficial to place observational constraints on this component. As mentioned before, this is a challenging task. Still, there are tantalizing hints of detections of WHIM gas. Wang & McCray (1993) found an excess emission component in ROSAT data around $`T2\times 10^6`$ K, which may be arising from this diffuse component. Soltan et al. (1996) detected an auto-correlation signature between the soft X-ray background and galaxies, as would be expected if warm-hot gas was distributed in large-scale structures. Even more exciting is a possible direct detection of diffuse soft X-ray emission associated with a filament of galaxies, by Scharf et al. (2000). Somewhat stronger evidence comes from a census of O VI absorbers at redshifts $`0.14z0.27`$ (Tripp, Savage & Jenkins (2000); Tripp & Savage (2000)), which together with conservative ionization corrections and metallicities implies $`\mathrm{\Omega }_{\mathrm{WHIM}}0.003_{0.002}^{+0.004}h_{75}^1`$, or $`\mathrm{\Omega }_{\mathrm{WHIM}}/\mathrm{\Omega }_b10`$% at the present epoch. However, some of these O VI absorbers may be photoionized, arising in cooler, low-density intergalactic gas at $`T10^4`$ K.

WHIM gas emission may be easiest to detect around high density regions such as clusters because that is where the density and temperature are highest within the WHIM range (cf. Figure 6), but this is not where the majority of WHIM gas is located. For instance, high-resolution Chandra spectra of regions between clusters that are free of bright sources may provide a detectable signal (Phillips et al., in preparation), and observations from XMM could possibly image filaments of warm-hot gas directly (Pierre, Bryan & Gastaud (2000)). Our simulations, however, predict the majority of WHIM gas is far away from galaxies and clusters, residing in the diffuse IGM. A promising avenue to detect this more typical WHIM gas is via absorption, as continuing observations with STIS aboard Hubble will detect many more O VI absorbers. Future X-ray satellites may be able to detect higher ionization absorbers such as O VII and O VIII that may also trace WHIM gas (Hellsten, Gnedin & Miralda-Escudé (1998)). The detection of this component is a key observational challenge, as the warm-hot intergalactic medium is rapidly become an integral part of our understanding of the evolution of baryons in the universe.

We thank Jeff Gardner, Ed Jenkins, Richard Mushotsky, Arielle Phillips, Jim Peebles and Todd Tripp for helpful discussions. RD is supported by NASA ATP grant NAG5-7066. RC and JPO are supported by NSF grants AST-9803137 and ASC-9740300. Support for GLB was provided by NASA through Hubble Fellowship grant HF-01104.01-98A from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. This work was supported by NASA Astrophysical Theory Grants NAG5-3922, NAG5-3820, and NAG5-3111, by NASA Long-Term Space Astrophysics Grant NAG5-3525, and by the NSF under grants ASC93-18185, ACI96-19019, and AST-9802568. Some of the simulations were performed at the San Diego Supercomputer Center. We also thank NCSA for use of their computing facilities.


format: text

text:
# Untitled Document

Interference of macroscopic superpositions

Italo Vecchi

Vicolo del Leoncorno 5 - 44100 Ferrara - Italy

email: vecchi@isthar.com

Abstract: A simple experimental procedure based on the Elitzur-Vaidman scheme is proposed to test the persistence of macroscopic superpositions. It is conjectured that its implementation will reveal the persistence of superpositions of macroscopic objects in the absence of a direct act of observation.

Schroedinger’s cats have kept physicists mulling since their appearance in $`[1]`$ and the properties of macroscopic superpositions are still hotly debated. The alleged impossiblity to detect superpositions of macroscopic objects (e.g. cats which are at the same time dead and alive) has given rise in recent years to a theory of environment-induced decoherence ( see $`[7]`$ for a detailed introduction to its different aspects and  for a survey of alternative approaches) which attributes the apparent absence of macroscopic superpositions to environment-induced perturbations. We propose here a simple experimental procedure based on Elitzur-Vaidman’s ”interaction-free measurement” to test the persistence of macroscopic superpositions. Our present experimental procedure yields different outcomes in the case when state vector reduction takes place (e.g. the cat is either dead or alive) and in the case where macroscopic superpositions persist. We conjecture that the implementation of our experiment will reveal the persistence of superpositions of macroscopic objects in the absence of a direct act of observation. Experiments on macroscopic quantum coherence have been proposed by other authors ( see $`[14]`$, $`[15]`$, $`[9]`$, $`[13]`$ ), but the argument proposed here appears to be new and quite feasible.

Interaction-free measurement is considered in the Elitzur-Vaidman scheme ($`[4]`$, see also $`[6]`$) and in earlier work by Dicke ($`[3]`$) and Renninger ($`[2]`$). In the Elitzur-Vaidman scheme, which has attracted considerable interest in recent years (see e.g. $`[8]`$), a robust interference technique is used to ascertain a system’s unknown macroscopic properties. The theoretical arguments motivating the present proposal are set forth in $`[5]`$, but they are not strictly necessary for an understanding of this paper. It may suffice to recall that in $`[5]`$ the assumptions underlying current decoherence theory are examined and any ”preferred” or ”pointer” basis is shown to depend on the observer and not to be an intrinsic property of the physical system.

Fig. 1

The formalism of decoherence theory encodes observer-dependent information into the system’s evolution prior to measurement. The fact that such encoding is spurious is demonstrated by the procedure proposed here, which implements a ”tilted basis” measurement along the lines proposed in the final remarks of $`[5]`$, so that interference terms are not averaged out by random phase changes (cf. $`[14]`$, par. 2). In this setting the Elitzur-Vaidman scheme reveals robust interference patterns that cannot be detected in the most obvious ”preferred basis”.

Let us consider now the device described in Fig. 1. The experimental setting, which is based on Mach-Zehnder interferometry, is introduced in $`[4]`$. The small circle in the lower left corner is a photon source, emitting one photon at a time. The arrows refer to the photon’s path while the thick hatched lines denote beam splitters, e.g. half-silvered mirrors. The dotted lines denote the possible positions of the absorber, as described below .

We adapt the procedure of $`[4]`$ by coupling the absorber $`F`$ to a quantum device $`QD`$, which creates a superposition $`\lambda _A|SA+\lambda _B|SB`$ of two states $`|SA`$ and $`|SB`$, corresponding to macroscopic events $`MEA`$ and $`MEB`$. In the original version of the Schroedinger’s Cat experiment $`QD`$ is a Geiger counter which is switched on for a given time and is connected to a relay that shatters a flask of hydrocyanidic acid. We may describe our quantum device $`QD`$ as a Geiger counter which is switched on for a given time and is connected to a relay that moves the absorber $`F`$. In Schroedinger’s original () argument one has $`\lambda _A^2=\lambda _B^2=\frac{1}{2}`$ and the macroscopic events $`MEA`$ and $`MEB`$ are the shattering of cyanide flask or its remaining intact, while $`|SA`$ and $`|SB`$ are the system’s states where the Cat is respectively dead or alive. In our setting too $`\lambda _A^2=\lambda _B^2=\frac{1}{2}`$, while $`MEB`$ corresponds to the insertion of the absorber $`F`$ on the photon’s path (Pos. B in Fig. 1) at time $`T=0`$, while $`MEA`$ corresponds to the absorber $`F`$ staying clear of the photon’s path ( Pos. A in Fig. 1). Accordingly $`|SA`$ corresponds to the absorber $`F`$ being in position $`A`$ while $`|SB`$ corresponds to $`F`$ being in position $`B`$.

We assume that during each experimental run the system $`S`$ beyond the wall $`W`$ is isolated in such a way that the only information accessible to the observer $`O`$ is channelled through the detector $`D`$. The system $`S`$ does not need to be fully physically isolated from the observer $`O`$ . Indeed the assumption that a system $`S`$ be completely isolated is impossible to realise in practice. What is needed in our setting however is just the absence of a coherent communication channel between $`O`$ and $`S`$. In other words in our setting during the experiment the observer-experimenter has no way to know what is going on in system $`S`$ behind the wall $`W`$, excepts for the readings of the detector $`D`$.

Both beam splitters $`BS1`$ and $`BS2`$ have transmission coefficient $`\frac{1}{2}`$. The state of a photon moving to the right is denoted by $`|1`$ and the state of a photon moving up is denoted by $`|2`$. We denote as $`|absorbed`$ the ”empty wave” corresponding to photon’s absorption. Reflection induces a change of phase of $`\frac{\pi }{2}`$ the photon’s wave function. The action of the beam splitter on the photon is ()

$$|1\frac{1}{\sqrt{2}}(|1+i|2),|2\frac{1}{\sqrt{2}}(|2+i|1)$$

while the action of the mirrors is

$$|1i|2,|2i|1.$$

We want to verify whether at time $`T=1`$, when the photon is fired, the system subsists as a superposition of a a state $`|SA`$ where the absorber $`F`$ is in position $`A`$ and of a state $`|SB`$ where the absorber in position $`B`$. There are two situations that we must consider. In the first situation the system has already collapsed into one of the ”macroscopically stable” states $`|SA`$ and $`|SB`$ so that one of the following two possibilities holds.

1a. The absorber $`F`$ is in position $`A`$. The process is described by

$$|1\frac{1}{\sqrt{2}}(|1+i|2)\frac{1}{\sqrt{2}}(i|2|1)$$

$$\frac{1}{2}(i|2|1)\frac{1}{2}(|1+i|2)=|1$$

so that the detector $`D`$ never clicks (cf. $`[4]`$).

1b. The absorber $`F`$ is in position $`B`$. The process is described by

$$|1\frac{1}{\sqrt{2}}(|1+i|2)\frac{1}{\sqrt{2}}(i|2+|absorbed)$$

$$\frac{1}{2}(i|2|1)+\frac{1}{\sqrt{2}}|absorbed$$

so that in this case the detector $`D`$ clicks with probability $`\frac{1}{4}`$.

Combining the above possiblities we see that, if the system $`S`$ has already collapsed into one of the ”macroscopically stable” states $`SA`$ and $`SB`$, then the probability that the detector $`D`$ clicks is

$$Pb(D)=\frac{1}{2}\frac{1}{4}=\frac{1}{8}.$$

This is the behaviour that we expect if state vector reduction has already taken place when the photon is fired, i.e. if the system’s state is a proper mixture of $`|SA`$ (case 1a) and $`|SB`$ (case 1b).

We consider now the case where the superposition of the states $`\frac{1}{\sqrt{2}}|SA`$ and $`\frac{1}{\sqrt{2}}|SB`$ corresponding to position $`A`$ and position $`B`$ subsists, but is perturbed by random interaction with the environment. As we will see, the measurement outcomes at the detector $`D`$ are different from those arising when either case 1a or case 1b holds.

2. The absorber $`F`$ is in a superposition $`\frac{a}{\sqrt{2}}|SA+\frac{b}{\sqrt{2}}|SB`$ of position $`A`$ and position $`B`$, where $`a`$ and $`b`$ are unknown uncontrollable phases such that $`|a|^2=|b|^2=1`$ . In this situation the photon crossing the superposed absorber $`F`$ undergoes a random phase change induced by the the coupling between the photon’s and the absorber’s states. The process is described by

$$|1\frac{1}{\sqrt{2}}(|1+i|2)$$

$$\frac{\alpha }{2\sqrt{2}}(i|2|1)+\frac{i}{\sqrt{2}}|absorbed\frac{1}{2}(|1+i|2)=$$

$$\frac{1}{2}\left[\left(\frac{\alpha }{\sqrt{2}}1\right)i|2\left(\frac{\alpha }{\sqrt{2}}+1\right)i|1\right]+\frac{i}{\sqrt{2}}|absorbed$$

where $`\alpha =\mathrm{cos}\theta +i\mathrm{sin}\theta `$ is the photon’s phase after it crosses the superposed absorber. Since the phase change is random, $`\theta `$ is spread uniformly on the interval $`[0,2\pi ]`$ . Thus the probability that the detector $`D`$ clicks is obtained averaging on the $`[0,2\pi ]`$ phase interval, i.e. on a randomly perturbed series of experimental runs

$$Pb(D)=\frac{1}{2\pi }_0^{2\pi }\left|\frac{\mathrm{cos}\theta +i\mathrm{sin}\theta \sqrt{2}}{2\sqrt{2}}\right|^2𝑑\theta =$$

$$\frac{1}{16\pi }_0^{2\pi }22\sqrt{2}\mathrm{cos}\theta +\mathrm{sin}^2\theta +\mathrm{cos}^2\theta d\theta =\frac{3}{8}$$

so that the measurement outcomes at the detector $`D`$ are different from those where the system behaves as if state vector reduction had already taken place. Heuristically we may describe the process in terms of ”photon fractions”, i.e. photon states with amplitude less than one, as follows. One half of the half photon following the upper path is absorbed by $`F`$ and the other half reaches $`BS2`$ where its $`|2`$ component combines with that of the half photon following the lower path. The standard (see $`[14]`$, par. 2) and plausible assumption is that the phases of the photon’s fractions which are not absorbed are randomly perturbed as the photon becomes entangled with the absorber. The result is just an instance of Feynman‘s rule that entangled particles do not exhibit interference other than with themselves ($`[11]`$, 1.6, III, see also $`[17]`$). In our case too ”the effect of photons being scattered is enough to smear out any interference effect” ($`[11]`$, 1.6-11). A random perturbation of the particle‘s phase blurs any interference pattern, so that superpositions appear to vanish. It should be stressed however that the assumption that the upper photon s phase is randomly parturbed as it crosses the absorber, while plausible and heuristically suggestive, is not an essential part of the argument. The key point here is that an entangled-perturbed particle does not exhibit interference with any other particle. It is the lack of interference between the upper scattered photon and the lower unperturbed photon that reveals the persistence of the absorber‘s superposition.

We observe that the persistence of macroscopic superpositions yields different outcomes also under the assumption that the photon’s phase is unaffected as it crosses the superposed absorber. Indeed, given a phase coefficient $`\alpha =\mathrm{cos}\theta +i\mathrm{sin}\theta `$, the probablity that the photon be detected at $`D`$ is

$$Pb(D)=\frac{1}{8}|32\sqrt{2}\mathrm{cos}\theta |\frac{32\sqrt{2}}{8}$$

where the equality holds only if $`\theta =0`$, i.e. if there is no coupling between the absorber’s and the photon’s phase. This case corresponds to a superposition of the solutions for 1a and 1b, but, since the initial conditions do not match those induced by the Geiger counter, the superposition principle is not applicable.

Finally we observe that, when the information on the sytem $`S`$ is extracted by the observer through the detector $`D`$, the absorber’s state no longer corresponds to a unitary evolution, so that thereafter the condition $`\lambda _A^2=\lambda _B^2=\frac{1}{2}`$ may be violated. The system must be re-inizialized after each measurement and then the procedure can be applied again. As long as the unitarity condition holds, i.e. as long as the device described here can be imbedded in a system whose evolution is unitary, the above argument holds. We refer to $`[10]`$, where a relevant argument pointing out the implications of the superposition principle at the macroscopic level is presented.

Our argument shows that in this setting entanglement prevents the system from duplicating the experimental outcomes of state-vector reduction. We conjecture that the experimental outcomes will confirm the persistence of macroscopic superpositions, as in case 2. A measured value $`Pb(D)\frac{1}{8}`$ would confirm the persistence of superpositions. An experimental outcome of $`Pb(D)=\frac{1}{8}`$ would provide strong evidence that state vector reduction has already taken place when the photon crosses the absorber, although it would not directly imply that this is the case, since, as shown above, there are values of $`\theta `$ which yield $`Pb(D)\frac{1}{8}`$, so that some peculiar distribution of phase changes corresponding to some coherent photon-absorber coupling might yield the same measurement outcome as state vector reduction.

It is worth stressing that, although the environmental states induced by the absorber are orthogonal, they do not correspond to the photon states which are being measured at the detector $`D`$. The measurement does not take place in the ”macroscopic” basis corresponding to the absorber’s position, but in the basis relative to the photon state at $`D`$, whereas the photon state which is being measured is just the projection of the system’s state vector on the photon’s state space at $`D`$. A heuristic example may clarify the issue. Suppose there are boys and girls in a room. The boys throw apples out of the window, the girls throw both apples and oranges. If the observer walks down the street and is hit by an orange he knows that it was a girl, if it is an apple he does not know whether it was a boy or a girl. Boys and girls are orthogonal, but the observer is measuring them through a tilted basis. In this case case absorber states are boys vs. girls, whereas photon states at $`D`$ are apples and oranges. We stress that the slanted correspondence between absorber and photon states is induced by the lower unperturbed photon. Indeed we can modify the device arranging the absorber so as to block the photon‘s upper path when it it is in Pos. A and the lower path when it is in Pos. B. Both the upper and the lower photon are scattered and become entangled and so they do not interfere. In this case the blurring of superposition does indeed mimic collapse, as in $`[14]`$.

One may consider the possibility that the photon phase is affected by the environment in such a way as to mimic collapse. The experimental outcomes for the device in Fig. 1 corresponding to collapse can be duplicated only assuming that the superposed absorber affects the lower photon’s phase in a very specific way. One has to assume that the lower photon is split into two equal amplitudes, one of which is left unperturbed, while the other undergoes exactly the same random changes that affect the phase of the upper photon as it crosses the absorber. It is easy to see that in this case the experimental outcomes are the same as in the case where state-vector reduction takes place: the perturbed photons cancel out and the unperturbed one is just what is needed to mimic collapse. Simply assuming that the superposed absorber induces a random perturbation of the lower photon’s phase is not enough to mimic collapse, since it yields the same outcomes as in our conjecture, i.e. with the upper and lower photon not interfering. The superposed absorber action on the far-away lower photon must be both strong and fine-tuned. The problem with this model is that it attributes to superpositions quite extraordinary physical properties, which distinguish them from ”ordinary” states. Since the difference between a superposition and an ”ordinary” (i.e. non-superposition) state depends only on the chosen basis this is quite perplexing, as the system would have to exhibit extraordinary physical properties in order to meet the observer’s expectations. In reality particle’s phases are quite sturdy and the can be changed only through strong interaction. It is precisely the relative insensitivity of particle’s phases to environmental perturbations that makes Zehnder’s interferometry possible.

It is worth pointing out that the argument presented in this letter appears related to Maris‘s recent work on fractional electrons in liquid helium $`[12]`$, which has attracted considerable interest. The phenomena examined in $`[12]`$ depend on a tilted-basis measurement, as superposed ”electrino bubbles” are generated by delocalised electrons (cf. $`[16]`$). The persistence of such bubbles, which carry a whole electron charge , is revealed in the photoconductivity experiments by ionic mobility measurements. The ionic mobilty of superposed ”electrino bubbles” is increased since they are smaller than non-delocalised bubbles. The second paragraph of $`[12]`$ deals with questions on quantum measurement theory that are related to those addressed here. The ”tilted-basis” measurement described here is characterised, similarly to those described in $`[12]`$ , by ”non-vanishing dispersion”, as pointed out in $`[16]`$.

References

 E. Schroedinger, Naturwissenschaften 23, 807 (1935).

 M. Renninger, Z. Phys. 158, 417 (1960).

 R.H. Dicke, Am. J. Phys. 49, 925 (1981).

 A. Elitzur and L. Vaidman, Found. Phys. 23, 987 (1993).

 I.Vecchi, at http://xxx.lanl.gov/abs/quant-ph/0002084 (see also: http://xxx.lanl.gov/abs/quant-ph/0102130).

 P.Kwiat, H. Weinfurter, T. Herzog and A. Zeilinger, Phys. Rev. Lett. 74(24), 4763 (1995).

 ”Decoherence and the Appearance of a Classical World in Quantum Theory” D.Giulini et al. ed. , Spinger, Berlin-Heidelberg-New York 1999.

 R. Penrose ”Shadows of the Mind”, Oxford University Press, Oxford, 1994.

 G.Ghirardi, at http://xxx.lanl.gov/abs/quant-ph/9810028 (1998).

 A.Bassi and G.Ghirardi, at http://xxx.lanl.gov/abs/quant-ph/0009020 (2000).

 ”The Feynman Lectures on Physics” , Addison-Wesley, 1963.

 H.J. Maris, Journal of Low Temperature Physics 120, 173.

 A.J. Leggett, Found. Phys. 29, 447 (1999)

 E.H. Walker et al. , Physica B, 339 (1988)

 F. Thaheld, Physics Letters A, 273, 232 (2000)

 R. Jackiw, C. Rebbi and J.R. Schrieffer, at http://xxx.lanl.gov/abs/cond-mat/0012370

 A. Zeilinger, Rev. Mod. Phys. 71, No.2 (1999).


text:
# Spinless impurities and Kondo-like behavior in strongly correlated electron systems

## I Introduction

Recently, magnetic properties induced by spinless impurities in correlated electron systems have attracted much interest. Especially, to probe antiferromagnetic correlations of High-$`T_c`$ cuprates the substitution of Cu sites with non-magnetic impurities such as Zn, Al, and Li, has been studied experimentally. According to NMR experiments, it was found that the substitution with spinless impurities induces local moments in the vicinity of impurities, which also show Kondo-like behaviors. For instance, the spin susceptibility in the vicinity of an impurity shows the temperature dependence like $`1/(T+T_\mathrm{K})`$, which implies the existence of the characteristic energy scale $`T_\mathrm{K}`$ analogous to Kondo temperature. Moreover the spin lattice relaxation rate $`1/T_1`$ shows Korringa-like behaviors, $`1/T_1TK^2`$, for $`T<T_\mathrm{K}`$. Here $`K`$ is the Knight shift. It is noted that in the vicinity of a spinless impurity the antiferromagnetic spin correlation which is developed in the bulk is much suppressed, and the magnetic correlation is dominated by the induced local moment. From theoretical points of view, it is non-trivial how this induced local moment governs the magnetic properties around an impurity, suppressing the antiferromagnetic correlation. In this paper, we shall deal with this issue. Although the experiments are carried out for High-$`T_c`$ cuprates which are essentially quasi-two-dimensional systems, it is expected that such effects may depend on the lattice structure and the dimensionality. Thus, we consider the Hubbard models with a spinless impurity in the spatial dimension $`D=1,2,`$ and $`3`$ to investigate how the dimensionality affects the induced magnetic properties. For $`D=1`$, the effects of a spinless impurity is incorporated into an open boundary condition as will be explained in the next section. Thus we consider the 1D Hubbard model with boundaries which is exactly solvable in terms of the Bethe ansatz method. We analyze the magnetic properties of this model using the exact solution and boundary conformal field theory. For $`D=2`$ and $`3`$, we derive the Korringa relation satisfied in the vicinity of a spinless impurity which is observed in NMR experiments. Our argument for $`D=2`$, $`3`$ is based upon the Fermi liquid theory in the presence of antiferromagnetic spin fluctuations.

The organization of this paper is as follows. In Sec. II, the 1D Hubbard model with a spinless impurity is considered. The spin susceptibility and the local density of states in the vicinity of an impurity are obtained based upon the Bethe ansatz exact solution. It is found that the induced moment is not screened completely at any finite temperatures. We also derive the spin-lattice relaxation rate $`1/T_1`$ which satisfies a relation analogous to the Korringa law. In Sec. III, we discuss about the 2D and 3D systems exploiting the antiferromagnetically correlated Fermi liquid theory. Summary is given in Sec. IV.

## II A spinless impurity in the one-dimensional Hubbard model

### A Mapping to the Hubbard model with boundaries and the Bethe ansatz exact solution

The effects of a single impurity in one-dimensional correlated systems have been extensively studied so far. If the interaction between fermions is repulsive, a potential scattering in the 1D Hubbard model is renormalized to an infinite strength, eventually, cutting the system into two half-infinite chains in the low-energy scaling limit. Thus at sufficiently low temperatures the system can be treated as the Hubbard chain with open boundaries, of which the hamiltonian is given by,

$$H=\underset{\sigma ,i=1}{\overset{L1}{}}c_{\sigma i}^{}c_{\sigma i+1}+h.c.+U\underset{i=1}{\overset{L}{}}n_in_i\mu \underset{\sigma ,i=1}{\overset{L}{}}n_{\sigma i}\frac{H}{2}\underset{i=1}{\overset{L}{}}(n_in_i)+V\underset{\sigma }{}n_{\sigma 1},$$
(1)
where the last term is a boundary potential. As we will see below, the low-energy spin dynamics around the impurity with which we are concerned are mainly described by this model, and the interaction or hopping between the two half-infinite chains is a subleading irrelevant interaction which can be incorporated by perturbative calculation.

The Bethe ansatz exact solutions of 1D correlated systems with boundaries have been studied by many authors. In connection with the spin dynamics in the vicinity of the boundary, an intriguing result was obtained for the supersymmetric t-J model by Essler. He obtained the divergent behavior of the boundary spin susceptibility as a function of a magnetic field $`H`$, i.e. $`\chi _{\mathrm{boundary}}1/H(\mathrm{ln}H)^2`$. It was first predicted by de Sa and Tsvelik that such a Curie-like behavior is universal for integral models with boundaries. Later, the similar behavior was also found for the Hubbard model at half-filling by Asakawa and Suzuki. In the next subsection, we shall show that this divergent behavior holds also for the case away from half-filling with finite $`U`$.

Here we summarize the basic equations which are relevant to the following arguments. The Bethe ansatz equations of the 1D Hubbard model with boundaries obtained by Schulz many years ago are,

$`e^{i2k_jL}e^{i\varphi _0(k_j)}`$ $`=`$ $`{\displaystyle \underset{\beta =1}{\overset{M}{}}}e_1(\mathrm{sin}k_j\lambda _\beta )e_1(\mathrm{sin}k_j+\lambda _\beta ),`$ (2)
$`{\displaystyle \underset{j=1}{\overset{N}{}}}e_1(\lambda _\alpha \mathrm{sin}k_j)e_1(\lambda _\alpha +\mathrm{sin}k_j)`$ $`=`$ $`{\displaystyle \underset{\genfrac{}{}{0pt}{}{\beta =1}{\beta \alpha }}{}}e_2(\lambda _\alpha \lambda _\beta )e_2(\lambda _\alpha +\lambda _\beta ),`$ (3)
where $`e_n(x)=\frac{x+\mathrm{i}nu}{x\mathrm{i}nu}`$, $`u=U/4`$, and $`\varphi _{0,L}`$ is a potential at boundaries. $`N`$ is the total number of electrons. $`M`$ is the total number of down spins. $`k_j`$ and $`\lambda _\alpha `$ are rapidities for charge and spin degrees of freedom, respectively. In the following, we consider only the case of repulsive boundary potentials. Thus the above equations have real roots. Putting $`k_j=k_j`$, $`\lambda _\alpha =\lambda _\alpha `$, and taking a continuum limit, we have the integral equations for the distribution functions of rapidities,

$`\rho (k)`$ $`=`$ $`{\displaystyle \frac{1}{\pi }}+{\displaystyle \frac{1}{\pi L}}\varphi _0^{}(k){\displaystyle \frac{1}{2\pi L}}{\displaystyle \frac{2u\mathrm{cos}k}{(\mathrm{sin}k)^2+u^2}}+\mathrm{cos}k{\displaystyle _B^B}{\displaystyle \frac{d\lambda }{\pi }}{\displaystyle \frac{u}{(\mathrm{sin}k\lambda )^2+u^2}}\sigma (\lambda ),`$ (4)
$`\sigma (\lambda )`$ $`=`$ $`{\displaystyle \frac{1}{\pi L}}{\displaystyle \frac{2u}{\lambda ^2+4u^2}}+{\displaystyle _Q^Q}{\displaystyle \frac{dk}{\pi }}{\displaystyle \frac{u}{(\lambda \mathrm{sin}k)^2+u^2}}\rho (k){\displaystyle _B^B}{\displaystyle \frac{d\lambda ^{}}{\pi }}{\displaystyle \frac{2u}{(\lambda \lambda ^{})^2+4u^2}}\sigma (\lambda ^{}).`$ (5)
$`N`$ and $`M`$ are given by,

$$_Q^Q\rho (k)𝑑k=\frac{2N+1}{L},$$
(6)
$$_B^B\sigma (\lambda )𝑑\lambda =\frac{2M+1}{L}.$$
(7)
Then the magnetization is expressed as,

$$\frac{S_z}{L}=\frac{1}{4}_Q^Q\rho (k)𝑑k\frac{1}{2}_B^B\sigma (\lambda )𝑑\lambda +\frac{1}{4L}.$$
(8)
The total energy is expressed in terms of the dressed energies,

$$\frac{E}{L}=_Q^Q𝑑k\left(\frac{1}{\pi }+\frac{1}{\pi L}\varphi _0^{}(k)\frac{1}{2\pi L}\frac{2u\mathrm{cos}k}{(\mathrm{sin}k)^2+u^2}\right)\epsilon _c(k)+_B^B\frac{d\lambda }{\pi L}\frac{2u}{\lambda ^2+4u^2}\epsilon _s(\lambda ),$$
(9)
where the dressed energies $`\epsilon _c(k)`$ and $`\epsilon _s(\lambda )`$ are determined by the integral equations,

$`\epsilon _c(k)`$ $`=`$ $`2\mathrm{cos}k{\displaystyle \frac{H}{2}}\mu +{\displaystyle _B^B}{\displaystyle \frac{d\lambda }{\pi }}{\displaystyle \frac{u}{(\mathrm{sin}k\lambda )^2+u^2}}\epsilon _s(\lambda ),`$ (10)
$`\epsilon _s(\lambda )`$ $`=`$ $`H+{\displaystyle _Q^Q}{\displaystyle \frac{dk}{\pi }}{\displaystyle \frac{u}{(\mathrm{sin}k\lambda )^2+u^2}}\epsilon _c(k){\displaystyle _B^B}{\displaystyle \frac{d\lambda ^{}}{\pi }}{\displaystyle \frac{2u}{(\lambda \lambda ^{})^2+4u^2}}\epsilon _s(\lambda ^{}).`$ (11)
If one fixes the magnetic field $`H`$, $`B`$ is determined by the equilibrium condition $`E/B=0`$, which is equivalent to the condition, $`\epsilon _s(B)=0`$. In the subsequent sections, we calculate the spin susceptibility and the local density of states using the above equations.

### B Spin susceptibility

In order to derive the spin susceptibility, we solve eqs.(4) and (5) for $`\sigma (\lambda )`$ using the Wiener-Hopf method, and obtain the magnetization, eq.(8).

Applying the Fourier transformation and shifting the argument, $`\lambda \lambda +B`$, we rewrite eq.(5) into,

$`\sigma (\lambda +B)=f_0(\lambda +B)+{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{d\lambda ^{}}{\pi }}R(\lambda \lambda ^{})\sigma (\lambda ^{}+B)+{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{d\lambda ^{}}{\pi }}R(\lambda +\lambda ^{}+2B)\sigma (\lambda ^{}+B),`$ (12)
where

$`f_0(\lambda +B)`$ $`=`$ $`{\displaystyle \frac{1}{L}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}{\displaystyle \frac{e^{u|\omega |}e^{i\omega (\lambda +B)}}{2\mathrm{cosh}u\omega }}+{\displaystyle _Q^Q}𝑑k{\displaystyle \frac{\rho (k)}{2\mathrm{cosh}\frac{\pi }{u}(\lambda +B\mathrm{sin}k)}},`$ (13)
$`R(x)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}{\displaystyle \frac{e^{u|\omega |}e^{i\omega x}}{2\mathrm{cosh}u\omega }}.`$ (14)
The last term of the right-hand side of eq.(12) is $`O(1/B^2)`$ for small magnetic fields. Thus we neglect it. Then eq.(12) can be solved by using the standard Wiener-Hopf method. The solution is expressed in terms of the following functions,

$`G^+(\omega )=\sqrt{2\pi }{\displaystyle \frac{(i\frac{u\omega }{\pi })^{i\frac{u\omega }{\pi }}}{\mathrm{\Gamma }(\frac{1}{2}i\frac{u\omega }{\pi })}}e^{\frac{iu\omega }{\pi }},`$ (15)
$`G^{}(\omega )=(G^+(\omega ))^1,`$ (16)
$`Q^+(\omega )+Q^{}(\omega )=G^{}(\omega )\stackrel{~}{f}_0(\omega ),`$ (17)
$`\stackrel{~}{f}_0(\omega )={\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑\lambda f_0(\lambda +B)e^{i\omega \lambda },`$ (18)
where $`Q^+(\omega )`$ ($`Q^{}(\omega )`$) is the analytic part of $`G^{}(\omega )\stackrel{~}{f}_0(\omega )`$ defined in the upper (lower) half plane. Fourier transforming eq.(12) and introducing the function $`\sigma ^+(\omega )=_0^{\mathrm{}}𝑑\lambda e^{i\omega \lambda }\sigma (\lambda +B)`$, we obtain the solution as, $`\sigma ^+(\omega )=G^+(\omega )Q^+(\omega )`$.

Now we derive $`Q^+(\omega )`$ as follows. For small magnetic fields, i.e. large $`B`$, and $`\lambda >0`$ the second term of $`f_0(\lambda +B)`$ is approximated as,

$$_Q^Q𝑑k\frac{\rho (k)}{2\mathrm{cosh}\frac{\pi }{u}(\lambda +B\mathrm{sin}k)}\frac{2N+1}{L}\frac{1}{2\mathrm{cosh}\frac{\pi }{u}(\lambda +B)}.$$
(19)
This driving term is essentially the same as the bulk contribution, with which we are not concerned. The first term of $`f_0(\lambda +B)`$ gives rise an interesting boundary effect. Using $`e^{u|\omega |}/2\mathrm{cosh}u\omega =_{n=1}^{\mathrm{}}(1)^{n1}e^{2nu|\omega |}`$ and the Laplace transformation,

$$\frac{2nu}{(\lambda +B)^2+(2nu)^2}=_0^{\mathrm{}}𝑑te^{(\lambda +B)t}\mathrm{sin}(2nut),$$
(20)
we rewrite the first term of $`f_0(\lambda +B)`$ as,

$$\frac{1}{L}_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{2\pi }\frac{e^{u|\omega |}e^{i\omega (\lambda +B)}}{2\mathrm{cosh}u\omega }=\frac{1}{\pi L}\underset{n=1}{\overset{\mathrm{}}{}}(1)^{n1}_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{2\pi }_0^{\mathrm{}}𝑑t\mathrm{sin}(2nut)\left[\frac{1}{\omega +it}\frac{1}{\omega it}\right]ie^{i\omega (\lambda +B)}.$$
(21)
The analytic property of eq.(21) solves eq.(17),

$$Q^+(\omega )=\frac{1}{L}\underset{n=1}{\overset{\mathrm{}}{}}_0^{\mathrm{}}𝑑t\mathrm{sin}(2nut)\frac{e^{tB}}{\omega +it}iG^{}(it)+\text{bulk terms}.$$
(22)
Finally, using eq.(8), we obtain the magnetization,

$$S_z=\frac{1}{2}_0^{\mathrm{}}𝑑\lambda \sigma (\lambda +B)=\frac{1}{2}\sigma ^+(0)\frac{1}{LB}+\text{bulk terms},$$
(23)
for large $`B`$, i.e. small magnetic fields. $`B`$ is related to $`H`$ from the condition $`\epsilon _s(B)=0`$. From eqs.(10) and (11), we have, $`H=Ce^{\pi B/2u}`$ for $`Hu`$. Here $`C`$ is an constant. Then the spin susceptibility $`\chi =S_z/H`$ behaves like,

$$\chi \frac{1}{L}\frac{1}{H(\mathrm{ln}H)^2}+\text{bulk terms}.$$
(24)
This $`H`$-dependence is the same as that found for the half-filling case. The above result implies that in 1D systems the magnetic moment induced by a non-magnetic impurity is not screened completely even at zero temperature. This behavior is analogous to the underscreening multi-channel Kondo effect, as pointed out by de Sa and Tsvelik. The leading $`H`$-dependence of eq.(24) is not altered, even if one includes irrelevant interactions such as the hopping between the two half-infinite chains.

In this section, we restrict our discussion to the zero temperature case. It is expected that at finite temperatures the boundary spin susceptibility behaves like $`\chi _{\mathrm{boundary}}1/T(\mathrm{ln}T)^2`$. In order to confirm this prediction, we need to explore thermodynamic Bethe ansatz method in the presence of boundaries. However, in the presence of boundaries, the entropy can not be expressed in terms of rapidity distribution functions in the continuum limit, because of the presence of spurious states for vanishing rapidities, and thus the usual technique of thermodynamic Bethe ansatz method is not applicable. If we limit the argument to sufficiently low temperature regions, undesirable contributions from the spurious state around the bottom of the energy spectrum may be small, and not give rise serious errors. Even if we admit this approximation, it is still a cumbersome task to solve thermodynamic Bethe ansatz equations numerically for low temperatures. Thus here we just give a field-theoretical argument to justify the above speculation. According to the boundary conformal field theory, the above divergent behavior of the spin susceptibility is due to the presence of a boundary entropy, $`S_{\mathrm{bound}}=T\mathrm{ln}(\sqrt{4\pi }R)`$. Here $`R`$ is the radius of the boson field of the Gaussian model which is the low-energy effective theory. If the leading irrelevant interaction is the marginal operator in the spin degrees of freedom, $`𝑱_L𝑱_R`$, we have $`RR_0g/\mathrm{ln}T`$ for small $`T`$. Then, the boundary spin susceptibility should behave like, $`\chi _{\mathrm{boundary}}1/T(\mathrm{ln}T)^2`$. Thus we expect that this temperature dependence which signifies the presence of an unquenched local moment may realize in this system.

### C Local density of states

In models solvable in terms of the Bethe ansatz method, the local density of states is defined as the derivative of the quantum number, which parameterizes rapidities, with respect to the pseudo-energy, i.e. $`I_j/\epsilon (k_j)`$. For the 1D Hubbard model, we can consider the density of states of holon and spinon, respectively. An interesting singular behavior due to the boundary appears in the spin degrees of freedom.

The local density of states of spinon as a function of energy is given by,

$$\rho _{\mathrm{spin}}(\epsilon )=\frac{\lambda }{\epsilon _s}\sigma (\lambda ).$$
(25)
In the absence of magnetic fields, $`B\mathrm{}`$, the solution of eq.(5) is expressed as,

$`\sigma (\lambda )={\displaystyle \frac{1}{L}}{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}{\displaystyle \frac{e^{u|\omega |}e^{i\omega \lambda }}{2\mathrm{cosh}u\omega }}+{\displaystyle _{\mathrm{}}^{\mathrm{}}}{\displaystyle \frac{d\omega }{2\pi }}{\displaystyle \frac{e^{i\omega \lambda }}{2\mathrm{cosh}u\omega }}{\displaystyle _Q^Q}𝑑k\rho (k)e^{i\omega \mathrm{sin}k}.`$ (26)
For $`\lambda 1`$, the first term of eq.(26) behaves like $`1/\lambda ^2`$, while the second term is just the order of $`O(e^{\pi \lambda /u})`$. Thus the main singular contribution comes from the former which is nothing but the boundary term. In a similar manner, from eq.(11) we obtain the asymptotic form of $`\epsilon _s(\lambda )`$ for large $`\lambda `$, i.e. $`\epsilon _s(\lambda )Ae^{\pi \lambda /u}`$, where $`A`$ is a constant. Then from eqs.(25) and (26), we have,

$$\rho _{\mathrm{spin}}(\epsilon )\frac{1}{\epsilon (\mathrm{ln}\epsilon )^2},$$
(27)
for small $`\epsilon `$. Thus the local density of states also shows the singular divergent behavior because of the presence of the boundary. It is noted that this result is similar to that of the underscreened multi-channel Kondo effect.

The important message of this and the previous subsections is that in 1D correlated electron systems the localized moment induced by a non-magnetic impurity is not quenched at any temperatures. The inclusion of irrelevant interactions such as hopping between semi-infinite chains does not change the result qualitatively. It should be stressed that this unquenched local moment is a particular property of the 1D systems where an impurity divides the system into two semi-infinite chains. Such a separation of the system is not possible in higher dimensional systems.

### D Spin lattice relaxation rate

Here we calculate the spin lattice relaxation rate $`1/T_1`$ in the vicinity of a spinless impurity, i.e. a boundary, using bosonization method and boundary conformal field theory. The same kinds of the calculations have been done for Heisenberg spin chains before. Some parts of the following results are similar to those obtained in ref.26 and 27. However combining them with the results from the Bethe ansatz exact solution, we shall see some new aspects. In the previous subsections, it was shown that the induced moment is not screened completely at any temperatures. Then one might expect that $`1/T_1`$ behaves like that of an isolated spin, $`1/T_1T\chi `$. However, as will be seen below, this naive expectation is incorrect.

According to the boundary conformal field theory, correlation functions for any operators in the vicinity of boundaries are obtained by the analytic continuation of the antiholomorphic part to the holomorphic part, $`O(z,\overline{z})=O_L(vt+ix)O_R(vtix)O_L(vt+ix)O_L(vtix)`$.  Following the standard technique, we have the asymptotic behaviors of the spin-spin correlation function in the presence of the boundary,

$`\chi (x,y,t)\left({\displaystyle \frac{\pi T}{v_s}}\right)^2[\mathrm{sinh}{\displaystyle \frac{\pi T}{v_s}}(xyvt)]^2`$ (28)
$`+e^{2ik_\mathrm{F}(xy)}{\displaystyle \underset{\nu =s,c}{}}\left[\left({\displaystyle \frac{\pi T}{v_\nu }}\right)^2{\displaystyle \frac{\mathrm{sinh}\frac{2\pi Tx}{v_\nu }\mathrm{sinh}\frac{2\pi Ty}{v_\nu }}{\mathrm{sinh}\frac{\pi T}{v_\nu }(x+y+v_\nu t)\mathrm{sinh}\frac{\pi T}{v_\nu }(x+yv_\nu t)\mathrm{sinh}\frac{\pi T}{v_\nu }(xy+v_\nu t)\mathrm{sinh}\frac{\pi T}{v_\nu }(xyv_\nu t)}}\right]^{\frac{K_\nu }{2}}.`$ (29)
Here $`v_s`$ and $`v_c`$ are the velocities of spinon and holon, respectively. $`K_c`$ is the Luttinger liquid parameter in the charge sector, and $`1/2K_c1`$. $`K_s=1`$ because of the $`SU(2)`$ symmetry of the spin sector. In the vicinity of the boundary, i.e. $`x,y,|xy|v_st`$, the staggered part (the second term) of eq.(29) is less relevant in comparison with the uniform part (the first term). Thus in contrast to the bulk behavior, the antiferromagnetic spin fluctuation is much suppressed, and the uniform part gives the dominant contribution to $`1/T_1`$ near the boundary. Fourier transforming eq.(29), we obtain, up to logarithmic corrections,

$$\frac{1}{T_1T}=\underset{\omega 0}{lim}\frac{1}{\omega }\underset{q,q^{}}{}\mathrm{Im}\chi (q,q^{},\omega )\frac{C}{v_s^2}+O(T^{2K_c}).$$
(30)
Here we have assumed that the hyperfine coupling constant is independent of $`q`$, and omitted it. $`C`$ is a temperature-independent constant. In the case that $`v_s`$ is a constant, the above result is equivalent to that obtained by Brunel et al. for Luttinger liquids with boundaries. The spinon velocity $`v_s`$ is related to the spin susceptibility obtained before, $`1/v_s=\chi _{\mathrm{bulk}}+\chi _{\mathrm{boundary}}/L`$. As claimed in the previous subsections, $`\chi _{\mathrm{boundary}}`$ should show enhanced local correlations like $`\chi _{\mathrm{boundary}}1/T(\mathrm{ln}T)^2`$. Thus near the boundary,

$$\frac{1}{T_1T}(\chi _{\mathrm{boundary}})^2.$$
(31)
Surprisingly, this relation is analogous to the Korringa relation. However it is noted that in contrast to the conventional Korringa law, the right-hand side of eq.(31) shows strong temperature dependence. As mentioned in the previous subsections, the induced local moment is not quenched completely. In spite of such an unscreened character of the moment, a la Korringa relation holds in the vicinity of a spinless impurity.

## III A spinless impurity in the 2D and 3D Hubbard models

In this section, we discuss the local magnetic properties caused by a spinless impurity in the 2D and 3D Hubbard model in the presence of bulk antiferromagnetic fluctuations, i.e. very close to the half-filling. The main purpose of this section is to derive the Korringa relation satisfied at the nearest neighbor of the impurity site, which is observed in the NMR experiment for cuprates. The model hamiltonian is given by,

$$H=\underset{k\sigma }{}E_kc_{\sigma k}^{}c_{\sigma k}+U\underset{i}{}n_in_i+E_0\underset{\sigma }{}n_{\sigma 0}.$$
(32)
Here $`E_k=2t_{a=1}^D\mathrm{cos}k_a`$ ($`D=2`$ or $`3`$), and the last term represents a spinless impurity localized at site $`0`$. We restrict our discussion to the case of the square lattice in 2D and the cubic lattice in 3D. Because of the presence of an impurity, correlation functions are non-local. In the case of $`U=0`$, the single particle Green’s function is given by,

$$G_{kk^{}}^0(\epsilon _n)=\frac{\delta _{kk^{}}}{i\epsilon _n+\mu E_k}+\frac{1}{i\epsilon _n+\mu E_k}\frac{E_0}{1E_0_{k^{\prime \prime }}\frac{1}{i\epsilon _n+\mu E_{k^{\prime \prime }}}}\frac{1}{i\epsilon _n+\mu E_k^{}},$$
(33)
where $`\mu `$ is a chemical potential. For $`U0`$, the single particle Green’s function is obtained by solving the equation,

$$\underset{k^{\prime \prime }}{}[(i\epsilon _n+\mu \epsilon _k)\delta _{kk^{\prime \prime }}\mathrm{\Sigma }_{kk^{\prime \prime }}(\epsilon _n)E_0]G_{k^{\prime \prime }k^{}}(\epsilon _n)=\delta _{kk^{}}.$$
(34)
The self-energy $`\mathrm{\Sigma }_{kk^{}}(\epsilon )`$ may be obtained by perturbative calculation in terms of $`U`$. However, in the following qualitative argument we do not need the explicit expression of $`G_{kk^{}}(\epsilon _n)`$.

Before discussing about magnetic properties, it is useful to sketch the spatial dependence of the density of states in the vicinity of an impurity. The density of states at the Fermi level is given by,

$$\rho (x,x^{})=(1/\pi )\underset{kk^{}}{}\mathrm{Im}G_{kk^{}}^R(0)e^{ikx}e^{ik^{}x^{}}.$$
(35)
To simplify the calculation, we consider the strong limit of an impurity potential, i.e. $`E_0U,t`$. The following argument do not change qualitatively even in the case of a finite $`E_0`$. Using eq.(33), we obtain the density of states at the Fermi level for the non-interacting system, $`U=0`$,

$`\rho _0(x,x^{})=N_{xx^{}}(\mu ){\displaystyle \frac{N_x(\mu )N_x^{}(\mu )}{N_0(\mu )}},`$ (36)
$$N_x(\epsilon )\underset{k}{}\delta (\epsilon E_k)e^{ikx}=_{\mathrm{}}^{\mathrm{}}\frac{ds}{2\pi }e^{is\epsilon }\underset{i=1}{\overset{D}{}}J_{n_i}(ts),$$
(37)
with $`J_n(x)`$, the Bessel function, and $`x=(n_1,n_2)`$ for $`D=2`$ and $`x=(n_1,n_2,n_3)`$ for $`D=3`$. Note that if $`x`$ or $`x^{}`$ is the impurity site, the density of states vanishes, $`\rho _0(x,0)=\rho _0(0,x^{})=0`$. We can easily show that if the electron density is close to the half-filling, i.e. $`|\mu |/t1`$, then $`N_x(\mu )O((\mu /t)^2)`$ for the site $`x`$ on the sublattice which includes the nearest neighbor site of the impurity, $`x_{\mathrm{n}.\mathrm{n}.}`$ (denoted by $`A`$-sublattice), and $`N_x(\mu )N_0(\mu )`$ for the site $`x`$ on the sublattice which includes the impurity site ($`B`$-sublattice). Thus from eq.(36) we immediately see that the local density of states around the impurity site shows strong spatial modulation similar to the Friedel oscillation. The period of the oscillation is $`1/k_\mathrm{F}`$ which is close to the half-filling value in this case. If $`x`$ and $`x^{}`$ belong to $`A`$-sublattice, the local density of states is nearly equal to that of bulk systems, $`\rho _0(x,x^{})\rho _0(|x|\mathrm{},|x^{}|\mathrm{})`$. On the other hand, if $`x0`$ and $`x^{}0`$ belong to $`B`$-sublattice, $`\rho _0(x,x^{})O((\mu /t)^2)`$. Since the density of states at the Fermi level is not renormalized by electron-electron interaction, this Friedel oscillation occurs for $`U0`$. This observation leads to an important implication for local magnetic properties on the nearest neighbor site of the impurity, $`x_{\mathrm{n}.\mathrm{n}.}`$. Because of the Friedel oscillation and the bipartite lattice structure, the local density of states on all the sites surrounding the site $`x_{\mathrm{n}.\mathrm{n}.}`$ is much suppressed provided that the impurity potential is sufficiently strong. Thus the spin on the site $`x_{\mathrm{n}.\mathrm{n}.}`$ is less screened than spins of electrons in the bulk. As a result, the local spin susceptibility on the site $`x_{\mathrm{n}.\mathrm{n}.}`$ is strongly enhanced.

Now we consider the spin lattice relaxation rate $`1/T_1`$ in the vicinity of a spinless impurity. For simplicity, we assume that the hyperfine coupling constant does not depend on $`q`$. We apply the general argument from the Fermi liquid theory to the case with a single spinless impurity. Then $`1/T_1`$ at the site $`x_i`$ is given by, up to constant factors,

$`{\displaystyle \frac{1}{T_1T}}`$ $`=`$ $`\underset{\omega 0}{lim}{\displaystyle \underset{q,q^{}}{}}{\displaystyle \frac{\mathrm{Im}\chi (q,q^{},\omega )}{\omega }}e^{iqx_i}e^{iq^{}x_i}`$ (38)
$`=`$ $`{\displaystyle \underset{q_1,q_2,q_3,q_4,k,k^{}}{}}\mathrm{Re}\mathrm{\Lambda }(q_1,k+q_2,k)\mathrm{Im}G_{kk^{}}^R(0)\mathrm{Im}G_{k+q_2k^{}+q_3}^R(0)\mathrm{Re}\mathrm{\Lambda }(q_4,k^{}+q_3,k^{})e^{iq_1x_i}e^{iq_4x_i},`$ (39)
where $`\mathrm{\Lambda }(q,k+q^{},k)`$ is a three point vertex function. The diagrammatic expression of eq.(39) which is the $`\omega `$-linear term of $`\mathrm{Im}\chi (q,q^{},\omega )`$ is shown in FIG.1. The detail derivation of this formula is given in ref.33 and 34. In the presence of strong antiferromagnetic fluctuations, it is plausible to assume that $`\mathrm{\Lambda }(q,k+q^{},k)`$ depends mainly on $`q`$ and $`q^{}`$. Thus in the following we discard the $`k`$-dependence of $`\mathrm{\Lambda }(q,k+q^{},k)`$. It is useful to rewrite eq.(39) in terms of quantities in the coordinate space,

$$\frac{1}{T_1T}=\pi ^2\underset{s,t}{}\rho (x_is,x_i+t)\rho (x_i+s,x_it)\mathrm{\Lambda }(x_i,x_i+s)\mathrm{\Lambda }(x_it,x_i).$$
(40)
Here,

$$\mathrm{\Lambda }(x,x^{})=\underset{q,q^{}}{}\mathrm{Re}\mathrm{\Lambda }(q,q^{})e^{iqx}e^{iq^{}x^{}}.$$
(41)
In the case that the site $`x_i`$ is far from the impurity, i.e. $`|x_i|a`$, where $`a`$ is a lattice constant, $`\rho (x_i+s,x_it)\rho (s+t)`$, $`\mathrm{\Lambda }(x_i,x_i+s)\mathrm{\Lambda }(s)`$, and thus the above expression is reduced to the usual formula of $`1/T_1`$ in bulk systems. It is also noted that if $`x_i=0`$, $`1/T_1T`$ vanishes, since $`\mathrm{\Lambda }(0,s)`$ includes $`G(0,x_j)`$ which vanishes as mentioned above. Here we are concerned with the case that the site $`x_i`$ is the nearest neighbor of the impurity site, $`x_i=x_{\mathrm{n}.\mathrm{n}.}`$.

To proceed further, we use a phenomenological expression for $`\mathrm{\Lambda }(q,q^{})`$. We assume that the three point vertex function $`\mathrm{\Lambda }(q,q^{})`$ consists of the part which is strongly enhanced by antiferromagnetic spin fluctuation and the local part which depends on $`q`$ and $`q^{}`$ weakly,

$$\mathrm{\Lambda }(q,q^{})\mathrm{\Lambda }_{\mathrm{AF}}(q)\delta _{q,q^{}}+\mathrm{\Lambda }_{\mathrm{loc}}(q,q^{}).$$
(42)
The antiferromanetically correlated part $`\mathrm{\Lambda }_{\mathrm{AF}}(q)`$ has a strong peak at $`q=Q`$, the staggered vector, and then approximated as, $`\mathrm{Re}\mathrm{\Lambda }_{\mathrm{AF}}(q)\mathrm{Re}\chi (qQ)=\chi (Q)/(1+(\xi _{\mathrm{AF}}(qQ))^2)`$. Here we used the phenomenological expression for $`\chi (qQ)`$. As mentioned above, $`\mathrm{\Lambda }_{\mathrm{loc}}(q,q^{})`$ is enhanced by local magnetic correlations at the site $`x_{\mathrm{n}.\mathrm{n}.}`$, i.e. $`\mathrm{\Lambda }_{\mathrm{loc}}(q,q^{})\mathrm{\Lambda }_{\mathrm{loc}}e^{iqx_{\mathrm{n}.\mathrm{n}.}}e^{iq^{}x_{\mathrm{n}.\mathrm{n}.}}`$ Then, eq.(39) is rewritten into,

$`\left({\displaystyle \frac{1}{T_1T}}\right)_{\mathrm{n}.\mathrm{n}.}`$ $``$ $`\left[{\displaystyle \frac{\chi (Q)}{(\xi _{\mathrm{AF}})^m}}\right]^2{\displaystyle \underset{k,k^{}}{}}\mathrm{Im}G_{k,k^{}}^R(0)\mathrm{Im}G_{k+Q,k^{}+Q}^R(0)`$ (44)
$`+[\mathrm{Re}\mathrm{\Lambda }_{\mathrm{loc}}]^2{\displaystyle \underset{k,k^{},q_2,q_3}{}}\mathrm{Im}G_{k,k^{},q_2,q_3}^R(0)\mathrm{Im}G_{k+q_2,k^{}+q_3}^R(0).`$
where $`m=2`$ for 2D systems and $`m=3`$ for 3D systems. Since $`\chi (Q)(\xi _{\mathrm{AF}})^2`$, the first term of eq.(44), which is the antiferromagnetically correlated part, is much suppressed compared to the second term, i.e. the local correlation part. Thus we obtain,

$$\left(\frac{1}{T_1T}\right)_{\mathrm{n}.\mathrm{n}.}(\mathrm{Re}\mathrm{\Lambda }_{\mathrm{loc}})^2.$$
(45)
Here we neglect all factors which are not enhanced by electron correlation. On the other hand, the local spin susceptibility at $`x_{\mathrm{n}.\mathrm{n}.}`$ is approximately given by, $`\chi _{\mathrm{loc}}\mathrm{Re}\mathrm{\Lambda }_{\mathrm{loc}}`$. Thus eq.(45) establishes the Korringa relation satisfied at the nearest-neighbor site of the impurity. As mentioned before, this relation is actually observed in NMR experiments.

## IV Summary and Discussion

We have discussed some magnetic properties analogous to the Kondo effect induced by a spinless impurity in strongly correlated electrons systems. In the 1D system, we have shown that the spin susceptibility and the local density of states near the impurity indicate divergent behaviors implying the presence of an unquenched local moment at any temperatures. We have also obtained the Korringa-like relation between the spin lattice relaxation rate and the local spin susceptibility, $`1/T_1T(\chi _{\mathrm{boundary}})^2`$. In the 2D and 3D systems, the antiferromagnetically correlated Fermi liquid theory has been applied. It has been shown that magnetic properties in the vicinity of a spinless impurity are dominated by the induced moment rather than the antiferromagnetic spin fluctuation developed in the bulk, and that the Korringa relation holds at the near neighbor site of the impurity.

The results obtained for the 1D system have an interesting implication to higher dimensional systems. Suppose a semi-infinite 2D Hubbard model with a boundary line, which is regarded as the coupled semi-infinite Hubbard chains. According to the results obtained in Sec. II, it is expected that at some finite temperatures the 1D-like strong spin correlations occurs in the vicinity of the boundary line leading to strongly enhanced density of states near the boundary. Such an enhanced electron correlation and one dimensionality of the boundary line may give rise strong fluctuations toward some surface phase transition. For instance, if there exists a pairing interaction in the bulk system, the paring correlation may be enhanced near the boundary leading to higher transition temperature than the bulk superconductivity. Actually, it is reported that $`\mathrm{Sr}_2\mathrm{RuO}_4`$ with lamellar microdomains of Ru metal shows the superconducting transition at the temperature higher than $`T_c`$ of the pure system, and that the superconductivity with higher-$`T_c`$ occurs in the vicinity of the boundary between $`\mathrm{Sr}_2\mathrm{RuO}_4`$ and Ru-metal. We would like to pursue this possible mechanism of the enhanced transition temperature in the near future.

###### Acknowledgements.

The author thanks N. Kawakami, H. Frahm, and K. Yamada for invaluable discussions. This work was partly supported by a Grant-in-Aid from the Ministry of Education, Science, Sports and Culture, Japan.