text:
# Ohmic Losses in Valence-band Photoemission Experiments

## I Introduction

Photoemission provides the most direct way to measure the the energy spectrum of electrons in a solid. The equation given in elementary modern physics texts in the section on the photoelectric effect is

$$E=E_\nu +ϵW,$$
(1)
where $`E`$ is the energy of the emitted electron, $`E_\nu `$ is the energy of the incident photon, $`ϵ`$ is the energy of the electron in the solid relative to the chemical potential and W is the work function. For monochromatic incident light, the spectrum of observed energies $`E`$ is the same as that for the electron energies $`ϵ`$ in the solid, except for a known shift.

This equation is known to be oversimplified in a number of respects. For example, an electron can scatter from other electrons on leaving the solid, giving up some of its energy, and even possibly producing secondary electrons which can escape from the sample. In fact, the mean free path of an electron moving with kinetic energies of 20 eV is short, of order $`10\AA `$, while the penetration depth of the light is much longer, meaning that only a small fraction of the electromagnetic energy absorbed actually leads to a photoemitted electron satisfying Eq. 1. Nonetheless, Eq. 1 is used for the analysis of the vast majority of angle-integrated photoemission spectroscopy (AIPES) experiments. The justification for this is the belief that the electrons that absorb a photon fall into two categories: (1) electrons that scatter strongly inelastically, losing an energy of several electron volts, and (2) electrons that do not scatter at all. Electrons in the first category either do not escape from the solid, or they do so with such a low kinetic energy that they can be discarded as uninteresting. Those in the second category are believed to follow Eq. 1. In view of the short mean free path, they must come from a thin layer near the surface, leading to the characterization of photoemission (PE) as a surface-sensitive probe.

Once stated in this way, a possible difficulty with this reasoning suggests itself. The argument depends on having a clean separation of energy scales. If an electron loses an intermediate amount of energy, the logic breaks down. In modern experiments, spectral features on a scale of 10 - 100 meV are very often of great interest. If many electrons have a substantial probability of undergoing inelastic scatttering with that much energy loss, the use of Eq. 1, which forms the basis of the interpretation of nearly all PE data, is no longer correct.

This paper will argue that this separation cannot always be be made, and that the assumptions underlying Eq. 1 can be dangerous in materials with relatively low conductivity or low effective carrier density. This includes many systems of current experimental interest.

Consider the photoemitted electron after it has left the solid at a speed $`v`$. Its electric field creates eddy currents and ohmic losses in the solid it is leaving behind. These losses are determined by the dielectric function $`ϵ(\omega )`$. If the spectral features are of width $`\mathrm{}\omega `$ 10 - 100 meV, then it is the Drude (or at least the low-frequency) dielectric response that determines the extent of the low-energy inelastic scattering. Dimensional considerations to be detailed below show that the probability for low energy scattering is

$$P_1=\frac{e^2}{v\mathrm{}}𝒜.$$
(2)
Here $`e^2/\mathrm{}v=(1/137)(c/v)0.8`$ for $`E_\nu =mv^2/2=20`$ eV. $`𝒜`$ is a dimensionless number that measures the absorption strength of the material in the Drude regime, but it is generally of order unity for isotropic materials with a Drude-like response. Thus the fraction of electrons that scatter is not small. The important issue is the typical energy loss of the electron. We will estimate below that a typical loss is of order $`\mathrm{}\sqrt{\sigma _0/\tau }`$ where $`\sigma _0`$ is the DC conductivity and $`1/\tau `$ is a relaxation rate - the width of the Drude portion of the function $`\sigma (\omega )`$. This is in the important 10-100 meV regime for materials in which the effective carrier density is low. This condition is fulfilled in many materials of current interest; in particular a material which undergoes a metal-insulator transition will usually pass through this regime.

Loss processes of the kind treated here, but generally at higher energy, have been considered before. For example, plasmon emission in X-ray PE causes satellite peaks in the measurement of core-level energies . Inelastic losses in EELS (electron energy loss spectroscopy) are due to essentially the same physical process . The theoretical description of these experiments is formally similar to that given here. The difference lies in the application of the concepts. The Drude part of the dielectric function is often temperature-dependent. This dependence can be dramatic in certain cases, as, for example, at a superconducting transition or a metal-insulator transition. This means that the extrinsic losses are also temperature-dependent. It is then not permissible to interpret temperature-dependent changes in an AIPES spectrum as being due only to changes in the density of states (DOS). The Drude part may also be sample-dependent. Thus, changes in AIPES from sample to sample may reflect changes in the conductivity, as well as changes in the DOS. Such effects are not so likely to occur in the plasmon part of the dielectric function. The low-frequency end of these theories has, however, never been of great interest before for PE. It is only the great advances in PE resolution, now $`20`$ meV for many experiments, which lends new interest to the application of these known concepts. A short summary of the results has appeared in previous publications .

In this paper we first give a classification of PE processes in order to situate the process treated here in the context of the overall theory of PE. Next, a detailed derivation of the loss function is presented, and its effect on observed spectra is illustrated in a simple case. Finally we give a short survey of some experimental systems that are candidates for the effect.

## II Classification of Loss Processes

The appropriate formalism is that due to Schaich and Ashcroft  and Caroli et al. . We present here a graphical interpretation of this formalism.

The conductivity (absorption of a low energy photon) in a solid is described by a current-current correlation. In the language of Feynman diagrams, we have a line representing an electron and a line representing a hole. Together, these make a loop with two vertices, the dipole matrix elements, at the ends. The PE experiment also involves the absorption of a photon, but there are two important differences between a conductivity diagram and a PE diagram for angle-resolved photoemission (ARPES). First, there is a measurement of the energy and momentum of the electron. Second, the electron line corresponds to a vacuum state, not a bulk state.

We represent the measurement by a rectangle inserted in the electron line. (In order to give formal expressions corresponding to each diagram, it is necessary to use the Keldysh rather than the Feynman formulation of perturbation theory, but the distinction does not affect the qualitative considerations given here.) This may now be called a current-current-current correlation function, as the electron line has been split into two. (The nomenclature is somewhat misleading since it would seem to imply that PE measures a third-rank tensor. It does not.)

The simplest diagram for a PE process is shown in Fig. 1. The electron propagates freely to the detector in this case. The measurement fixes the energy and momentum of the electron line, so it contributes only a delta function to the diagram and what we have is the textbook process of Eq. 1. If the momentum is averaged over, then the experiment is AIPES, which is the experiment we will concentrate on in this paper.

Fig. 2 is the most-studied of all the diagrams. The thick solid line is the full propagator for the hole. This diagram represents intrinsic many-body effects. All many-body effects on the propagation of a hole in the bulk system must be included, in principle. The reason for calling this diagram ”intrinsic”, is that nothing about the photoemission process itself is involved. The emitted electron still propagates freely to the detector. The diagram of Fig. 2 includes that of Fig. 1. To the extent that only these diagrams are important, PE measures the spectral function of a hole. In simpler language, it measures the occupied DOS.

Fig. 3 is the effect of the interaction of the emitted electron with the hole left behind. This is usually called a final state interaction, an ”extrinsic” process, as it depends on the PE process itself. The rate of such interactions decreases with incoming photon energy, and this diagram is usually neglected. This is called the sudden approximation. Fig. 4 represents an extrinsic process in the Born approximation. The thick interaction line is fully screened, and thus includes interactions with all the electrons in the material. This changes the energy and the momentum of the observed particle. This diagram is the subject of the present work. Note that the distinction between the final state interactions and the extrinsic inelastic loss processes is a clean one. The two should not be confused. The ohmic losses are, roughly speaking, due to interactions with the dynamic image charge, not the photohole.

This is not a complete listing of all possible processes, as one can also have scattering from impurities of the photoelectron, reflection at the surface, and so on. We have tried to give a reasonable overview of the possibilities that make an interesting comparison or contrast with the ohmic process.

## III Derivation for Cubic Systems

We begin with the general problem of a sample that occupies the half-space $`z<0`$. There is a time-dependent charge density $`\rho (\stackrel{}{r},t)`$ outside it. The Fourier component at frequency $`\omega `$, $`\rho (\stackrel{}{r},\omega )`$, sets up an electric field for $`z<0`$ given by

$$\stackrel{}{E}(\stackrel{}{r},\omega )=_{z^{}>0}\frac{\rho (\stackrel{}{r}^{},\omega )}{|\stackrel{}{r}\stackrel{}{r}^{}|}\frac{d^3r^{}}{1+ϵ(\omega )}.$$
(3)
This result follows from the usual image-charge calculation, namely the solution of Maxwell’s equations with the standard surface boundary conditions. Nevertheless, there are certain approximations. The normal skin effect is assumed because the wavevector dependence of $`ϵ`$ is neglected. This could break down for clean materials or at quite high frequencies. Fortunately, both of these conditions are outside of the applications which interest us here. Secondly, the real surface is not of infinitesimal thickness, and corrections for this have been neglected. Finally, we have specialized to the case that the material has cubic symmetry: the electric field and the electric displacement at any given frequency are simply proportional. The field drives currents that generate heat. The total energy thus produced is:

$$Q=\frac{1}{2}_{z<0}𝑑\omega \mathrm{}\sigma (\omega )|\stackrel{}{E}(\stackrel{}{r},\omega )|^2d^3r.$$
(4)
These results are quite general for any charge density for which the integral in Eq. 4 is convergent. Now we specialize to the PE process, which is modeled as a charge generated at the surface at time $`t=0`$ leaving the solid at an angle $`\theta `$ to the normal. The charge density is

$$\rho (\stackrel{}{r},t)=e\mathrm{\Theta }(t)\delta (zvt\mathrm{cos}\theta )\delta (y)\delta (xvt\mathrm{sin}\theta ),$$
(5)
where $`\mathrm{\Theta }(t)`$ is the step function. Since in most experiments $`\theta \pi /10`$, we set $`\mathrm{cos}\theta =1`$ henceforth. The Fourier transform is

$$\rho (\stackrel{}{r},\omega )=\frac{e}{2\pi v\mathrm{cos}\theta }\delta (x)\delta (y)e^{i\omega z/v}\mathrm{\Theta }(z).$$
(6)
Using Eq. 3, we have

$$\stackrel{}{E}(\stackrel{}{r},\omega )=\frac{e}{2\pi v\mathrm{cos}\theta }\frac{2}{1+ϵ(\omega )}_0^{\mathrm{}}𝑑z^{}e^{i\omega z^{}/v\mathrm{cos}\theta }\frac{\stackrel{}{r}z^{}\widehat{z}}{|\stackrel{}{r}z^{}\widehat{z}|^3}$$
(7)
as the field inside the sample. Substituting this into Eq. 4 yields

$$Q=\frac{2e^2}{\pi v}𝒞_0^{\mathrm{}}𝑑\omega \frac{\mathrm{}\sigma (\omega )}{\omega |1+ϵ(\omega )|^2},$$
(8)
where $`𝒞`$ is a dimensionless integral:

$$𝒞=_0^{\mathrm{}}\rho 𝑑\rho _0^{\mathrm{}}𝑑z\left(\rho ^2|I_1(\rho ,z)|^2+|I_2(\rho ,z)|^2\right),$$
(9)
with

$$I_1=_0^{\mathrm{}}\frac{e^{i\zeta }d\zeta }{\left[\rho ^2+(z+\zeta )^2\right]^{3/2}}$$
(10)
and

$$I_2=_0^{\mathrm{}}\frac{e^{i\zeta }(z+\zeta )d\zeta }{\left[\rho ^2+(z+\zeta )^2\right]^{3/2}}.$$
(11)
$`𝒞`$ may be estimated as about 2.6.

This classical calculation yields not only the total energy loss but also the loss into each frequency interval $`d\omega `$. Quantum-mechanically, this is interpreted as the inelastic energy loss at this energy due to a scattering event from the bulk electrons. It was calculated assuming no change in the trajectory of the electron, and is therefore essentially the Born approximation. Fig. 4 could also have several scattering lines. This would correspond to approximations beyond Born. When such processes are important, then the amount information contained in a PE spectrum decreases, just as in the case of multiphonon events in neutron scattering (for example).

The relative differential probability to lose energy $`\mathrm{}\omega `$ is

$$P(\omega )=\frac{2e^2𝒞}{\pi \mathrm{}v\omega ^2}\frac{\mathrm{}\sigma (\omega )}{|1+ϵ(\omega )|^2}.$$
(12)
There is also a probability for forward scattering $`P_0`$. This is the probability that an electron loses zero energy. Including this possibility, the observed intensity is given by

$$I(\omega ,T)=P_0(T)N(\omega )f(\omega )+_{\mathrm{}}^{\mathrm{}}P(\omega \omega ^{},T)N(\omega ^{})f(\omega ^{})𝑑\omega ^{},$$
(13)
Here $`N(\omega )`$ is the DOS and $`f`$ is the Fermi function. The total normalization is given by our assumption that there is either one or zero scatterings:

$$1=P_0+_0^{\mathrm{}}P(\omega )𝑑\omega .$$
(14)
Equations having the form of Eq. 13 are well-known in considerations of background in PE. For example, the Shirley background would result from the phenomenological assumption that $`P(\omega )`$ is a constant to be fit to experiment .

## IV Experimental consequences

¿From the point of view of experiment, one would like rules of thumb for recognizing the situations in which the ohmic effect is likely to be important. To this end, we investigate more closely the function $`P(\omega )`$. Note first that

$$P(\omega )=\frac{2e^2𝒞}{\pi \mathrm{}v}\frac{\mathrm{}\sigma (\omega )}{\omega ^2|1+ϵ(\omega )|^2}$$
(15)
vanishes at small frequencies if the system is a true insulator, since $`\mathrm{}\sigma `$ vanishes faster than any power of $`\omega `$. This must be true from a physical standpoint: an insulator cannot absorb electromagnetic energy. If the system is a very good metal, then $`\sigma _0\mathrm{}\sigma (\omega =0)`$ is large, and at low frequencies $`ϵ(\omega )4\pi i\sigma _0/\omega `$ leading to

$$P(\omega )\frac{1}{16\pi ^2\sigma _0}\frac{2e^2𝒞}{\pi \mathrm{}v}.$$
(16)
The vanishing of $`P(\omega )`$ in this limit is also required by elementary physics. A perfect metal must completely screen out the field at sufficiently low frequencies. We conclude that materials in the intermediate range are of most interest. Certainly if PE is done across a metal-insulator transition this range is explored.

There is no general formula that determines the importance of ohmic losses, because $`ϵ(\omega )=1+4\pi i\sigma /\omega `$ can take many forms, with the consequence that even materials with similar DC conductivities may give quite different results. Furthermore, the effect depends on the frequency range which is being explored. However, it very often happens that in the low-frequency regime the two-parameter Drude function

$$\sigma (\omega )=\frac{\sigma _0}{1i\omega \tau }$$
(17)
is a reasonable approximation. In this paper, we will focus on this particular form. Even for materials for which the Drude form is not a very good quantitative approximation for the conductivity, it may still be true that one parameter for the overall scale and one for the width are sufficient. In this case, the results below should give a reasonable guide.

The Drude form for $`P(\omega )`$, obtained by substituting Eq. 17 into Eq. 12 and rearranging, is:

$$P(\omega )=\frac{2e^2𝒞}{\pi \mathrm{}v}\frac{\sigma _0\tau ^2}{4}\left\{\left[\omega ^2\tau ^2(2\pi \sigma \tau 1/2)\right]^2+2\pi \sigma _0\tau 1/4\right\}^1.$$
(18)
The first thing to note about this function is its normalization. Integrating Eq. 18 we find

$$_0^{\mathrm{}}P(\omega )𝑑\omega =\frac{e^2𝒞}{4\pi \mathrm{}v}.$$
(19)
In terms of the notation of Eq. 2, this says $`𝒜=𝒞/40.65`$. Substituting in a velocity corresponding to the energy $`E_\nu =20`$ eV we find $`P_1=0.5`$. Hence the integrated effect is substantial. Since this result is independent of the parameters, it is a generic result not limited to the particular details of the Drude form. Bear in mind that this is not the only contribution to $`P_1`$ since losses at higher energies, for example losses due to interband processes, will also contribute. On the other hand, the fact that $`P_1`$ is comparable to unity at typical experimental energies also implies that we near the limit of validity of the Born approximation.

Written in terms of the variable $`\omega ^2\tau ^2`$, $`P(\omega )`$ is a Lorentzian centered at $`2\pi \sigma \tau 1/2`$ and width $`\sqrt{2\pi \sigma _0\tau 1/4}`$. In terms of the physical variable $`\omega `$, we have two regimes. If $`2\pi \sigma _0\tau <1/2`$, then the function is peaked at $`\omega =0`$, and falls off at large frequencies as $`\omega ^4`$ on a characteristic scale of $`1/\tau `$. If $`2\pi \sigma _0\tau >1/2`$, then the function is peaked at

$$\omega _{peak}=\sqrt{2\pi \sigma _0/\tau 1/2\tau ^2}\sqrt{2\pi \sigma _0/\tau },$$
(20)
and the width of the peak is

$$\mathrm{\Delta }\omega =\frac{1}{\tau }\left(2\pi \sigma _0\tau 1/4\right)^{1/4}\left(\frac{2\pi \sigma _0}{\tau ^3}\right)^{1/4}.$$
(21)
The limiting behaviors for $`2\pi \sigma _0\tau >>1`$ are also given.

The overall behavior of $`P(\omega )`$ is illustrated in Fig. 5. It is plotted for a fixed value of $`\sigma _0/\tau `$. In the Drude model

$$\frac{\sigma _0}{\tau }=\frac{ne^2}{m^{}}=\frac{\omega _p^2}{4\pi },$$
(22)
where $`\omega _p`$ is an effective plasma frequency. Thus, in a system where the effective carrier concentration is fixed and the temperature dependence comes either from the relaxation time or from $`\sigma _0`$, Fig. 5 represents the evolution of $`P(\omega )`$ with temperature. This graph is for a low effective carrier density: $`\sigma _0/\tau =10^4`$ meV2, which corresponds to an effective plasma frequency of 354 meV. The relaxation times are in units of meV-1 where $`\mathrm{}=1`$, that is $`6.6\times 10^{13}s`$. As $`\tau `$ increases, the peak position moves from $`\omega =0`$ to the limiting value $`\sqrt{2\pi \sigma _0/\tau }`$, while it continually sharpens up.

These structures in $`P(\omega )`$ manifest them selves in the observed intensity according to Eq. 13. Consider the case of zero temperature. Then $`f(\omega )`$ becomes $`\theta (\omega )`$ and we have:

$$I(\omega )=N(\omega )\theta (\omega )P_0+_\omega ^0P(\omega ^{}\omega )N(\omega ^{})𝑑\omega ^{}.$$
(23)
We are only interested in the second term, the inelastic component, here. Furthermore, we shall take a constant DOS $`N(\omega )=N_0=1`$ so that the features in the observed intensity come only from the energy dependence of the scattering. Thus we normalize $`P(\omega )`$ to unity.

$$I_{inelastic}(\omega )=_\omega ^0P(\omega ^{}\omega )𝑑\omega ^{}$$
(24)
The integral rearranges the apparent DOS. The DOS will appear unchanged at frequency $`\omega `$ if $`|\omega |>>\omega _{peak}`$ where $`\omega _{peak}=0`$ or $`\sqrt{\frac{2\pi \sigma _0}{\tau }\frac{1}{2\tau ^2}}`$ depending on the size of $`\sigma _0\tau `$. The photoelectrons observed at large $`|\omega |`$ are shifted, but it appears as if the (flat) DOS is unaffected at these frequencies. At low frequencies the situation is more complicated. The basic scenario for a constant $`\sigma _0/\tau `$ (effective carrier density) is this: at ’small’ $`\tau <<\tau _c`$ $`P(\omega )`$ is peaked around $`\omega =0`$ so the PE loss is negligible. At the crossover lifetime $`\tau _c=\sqrt{\frac{m^{}}{4\pi ne^2}}`$ the peak of $`P(\omega )`$ begins to shift to a higher frequency. In addition, $`P(\omega )`$ is very spread out for $`\tau `$ in this region, so the apparent loss becomes large and may appear pseudogap. When $`\tau >>\tau _c`$ the shape of the observed intensity $`I(\omega )`$ looks unchanged, but the entire curve is shifted downwards by an amount comparable to the effective plasma frequency, and this may look like a gap. This evolution of $`I(\omega )`$ with increasing $`\tau `$ is shown in Fig. 6. Bear in mind that this is only the inelastic component and it must be combined with the elastic intrinsic part. This may involve fitting since the overall normalization of the elastic part $`(P_0)`$ would normally be hard to determine.

It should be noted that while most of the peak shift of $`P(\omega )`$ occurs over a relatively small range of $`\tau `$ the evolution of $`I(\omega )`$ takes place over several orders of magnitude of $`\tau `$ since $`P(\omega )`$ is spread out over a large frequency range. We can see this by plotting the shift of the intensity curve vs. $`\tau `$ as in Fig. 7. The shift is defined by the position in energy of the point where the intensity is half its maximum, a popular definition of the gap in PE experiments.

## V Survey of Experimental Systems

Experimentally, the best candidates for gap-like structures are the one-dimensional conductors. There are several one-dimensional materials that show the development of some sort of a gap in the spectrum without corresponding observations in other properties. We mention here only two.

TTF-TCNQ is a one-dimensional organic chain compound. It has a $`2k_F`$ charge-density-wave transition at $`T_P=54K`$ with an energy gap $`2\mathrm{\Delta }40`$ meV, as measured by the activation energy of transport and thermodynamic quantities. The photoemission spectrum  is very much at odds with the transport and thermodynamics. At $`70K`$, there is a gap of about $`0.150.2eV`$, judging by the peak in an energy scan at an angle that corresponds to $`k_F`$. This peak does not move toward the Fermi energy as T increases. Rather, the spectrum just broadens out, but without increasing the weight at the Fermi energy. In fact, the weight if anything decreases . This material is a good candidate for a pseudogap of partly extrinsic origin.

(TaSe4)2I presents a somewhat similar scenario. It has a $`2k_F`$ charge density wave transition at about $`250K`$. A number of its properties: resistivity, magnetic susceptibility, and optical conductivity can be explained by assuming that the gap is roughly 250 meV, and remains at about this value, independent of temperature, even after long-range positional order is lost . But the gap in photoemission is about 500 meV, independent of temperature . Here we may have a gap which is split about 50-50 between intrinsic and extrinsic.

The theory presented above for a Drude-like conductivity is not applicable to the one-dimensional case as there is certainly a genuine pseudogap or gap, depending on temperature, and the systems are not cubic. The factor of two between the PE gap and the gap inferred from other properties in (TaSe4)2 We would find a natural explanation if the extrinsic loss turns on at the intrinsic gap frequency.

The colossal magnetoresistance manganate materials are good candidates for extrinsic losses in the neighborhood of their metal-insulator transition. The cubic material La0.67Ca0.33MnO3 is in fact the only system for which model calculations of the extrinsic losses have been carried out . The material has a metal-insulator transition at $`260K`$. In the ’metallic’ state at $`80K`$, there is a strong negative slope in $`I(\omega )`$ for at least 0.6 eV below $`\mu `$ . There is a sharp break in slope at $`\mu `$, presumably indicative of a nonzero density of states at $`\mu `$. But in the insulating state at $`280K`$ there appears to be no Fermi edge at all - the observed intensity is flat at $`\mu `$ and weight has moved back from $`\mu `$. These features might be taken as indicating the presence of a pseudogap which which opens in the insulating state. However, this movement of spectral weight can be produced by extrinsic effects, as shown in .

The layered material La1.2Sr1.8Mn2O7, on the other hand, has been shown by Dessau and Saitoh  to have a momentum-dependent pseudogap somewhat similar to the high-Tc materials, but rather larger, in the range of $`200meV`$ or so at the maximum. Again, the momentum dependence indicates that the pseudogap is mostly of intrinsic origin. The layered manganate mterials are an interesting case because of their extremely high resistivities. La1.2Sr1.8Mn2O7 itself has an a-b plane resistivity of about 6 m $`\mathrm{\Omega }`$-cm at T = 50 K, which is in the ’metallic’ phase. This is an extraordinarily high value. It is so close to being insulating that one would probably not expect ohmic effects .

Another class of poorly conducting materals are the Kondo insulators. This class of systems shows a resistivity dominated by the Kondo effect at higher temeratures crossing over to activated behavior at low temperatures, with absolute values which can be quite high. The resistivity as a function of temperature is monotonic, in contrast to the CMR systems. For example, the system YbB12 has such a crossover at about $`50K`$. This material also shows some indications of extrinsic processes. Even on the metallic side, the observed density of states slopes downward strongly to the Fermi energy . As the resistivity increases, there is a further depression of weight near $`E_F`$, with finally a gap-like structure of about 20 meV appearing at 6 K.

The best-known pseudogap of all is that in the high-temperature superconductor Bi2Sr2CaCu2O8 . In this case, the pseudogap is both momentum- and temperature-dependent, appearing mainly along the Cu-O bond direction. This is a strong indication that the pseudogap is largely intrinsic, as the extrinsic losses should not have a strong momentum dependence. Bi2Sr2CaCu2O8 is roughly tetragonal and the theory developed for the cubic case is not applicable. Its conductivity in the a-b plane is fairly high. The c-axis conductivity, which is very low, also plays a role, however. The theory for such a strongly anisotropic system has not been worked out in detail, and the combined effect of $`ϵ_{xx}`$ and $`ϵ_{zz}`$ is not yet completely clear. it seems most likely that lineshapes are affected by extrinsic effects, but that the size of the pseudogap is correctly given by the naive analysis.

The high-temperature superconductors and related cuprates are a class of systems for which some energy-dependent PE has been done. This is very interesting from the point of ohmic effects, which are predicted to decline as $`1/v1/\sqrt{E_\nu }`$, as seen from Eq. 12. For example, Sr2CuO2Cl2 has been investigated using ARPES for $`E_\nu =22eV`$ and $`E_\nu =35eV`$  and Pb-doped Bi2Sr2CaCu2O8 at $`E_\nu =32eV`$, $`E_\nu =40eV`$ and $`E_\nu =50eV`$ . In every case, changes in intensity as a function of $`E_\nu `$ are seen. In the case of Pb-doped Bi2Sr2CaCu2O8 they go in the direction of increasing binding energy as $`E_\nu `$ is decreased, which is in qualitative agreement with the present scenario. In addition, the lineshapes are rather broad and also energy dependent.

It is evident from this brief survey that an important direction for development of the theory would be a generalization to tetragonal and orthorhombic materials, as relatively few of the most interesting systems are cubic. The image calculation can be performed analytically for the tetragonal case if the emission is along the c-axis, but the loss function is rather complicated and will be published elsewhere. The orthorhombic case is even more complicated, as the image method is not applicable.

## VI Conclusion

Changes in the density of states of electronic systems on the scale of 10 to 100 meV are of great interest in today’s PE experiments. These energies often coincide with the scale of relaxation rates and conductivities or combinations thereof. This means that ohmic losses may affect the interpretation of the data. This effect is distinct from others, such as final state effects, that have been considered in detail in the past.

General formulas for the loss in terms of the dielectric function may be given. In the case where a Drude form is appropriate, we may characterize the situation as follows. The inelastic loss will be most pronounced if the lifetime is roughly $`\tau _c1/(2\pi \sigma _0)`$. The magnitude of the loss is on the scale of the maximum peak shift $`\omega _{peak}=\sqrt{2\pi \sigma _0/\tau }`$ near the effective plasma frequency. If the effective plasma frequency is large, as in a good metal, then the inelastically scattered electrons are emitted well below the chemical potential as secondaries, and the elastic peak is unsullied. If the product $`\sigma _0\tau `$ is very small, as in a near-insulator, the scattered electrons are emitted with very small loss and cannot be distinguished from the unscattered ones. In the intermediate case, however, a combination of the elastic and inelastic terms is required to give the observed spectrum. This ohmic effect may account for the observation of a gap or pseudogap in some systems with relatively low intrinsic conductivities in one or more directions. In general, however, low conductivity is a necessary, not a sufficient, condition, for the ohmic effects to be important. The effective carrier density and the frequency range of the measurements must also be taken into account.

## VII Acknowledgements

We would like to thank M. Rzchowski, M.B. Webb, A.V. Chubukov, L. Bruch, M. Norman, J.‘Joyce, and F.J. Himpsel for helpful discussions. This work was supported by the NSF Materials Theory Program, Grant No. DMR-0081039, and by the Materials Research Science and Engineering Center Program, Grant No. DMR-96-32527.


format: text

text:
# Orientation field model for chiral branching growth of bacterial colonies

## Abstract

We present a new reaction-diffusion model for chiral branching growth of colonies of the bacteria Paenibacillus dendritiformis. In our model the bacteria are represented by a density field with non-linear diffusion and a complex scalar field which represents bacterial orientation. The orientation field introduces anisotropy into the flux of bacteria, representing self-propulsion along their long axes. The model can also reproduce tip-splitting growth of other strains (shorter bacteria) of the same species. The model can capture changes of small number of bacteria, thus it can be used to study the open question of transitions between tip-splitting and chiral dendritic growth.

Bacteria display various chiral properties. Mendelson et al.  showed that long cells of Bacillus subtilis can grow in helices, in which the cells form long strings that twist around each other. They have also showed that the chiral characteristics affect the structure of the colony. Matsuyama et al.  have found that colonies of B. subtilis can grow into tip-splitting patterns with a global rotation about the center of the colonies (global twist). Ben-Jacob et al.  have found that some strains of Paenibacillus dendritiformis (see Ref.  for identification of the bacteria) exhibit similar patterns with global twist. They have also found that other strains (referred to as $`𝒞`$ (chiral) morphotype) present different chiral property, chiral dendritic growth with local twist (which they refer to as strong chirality).

Colonial patterns of $`𝒞`$ bacteria grown on semi-solid agar are characterized by chiral dendritic branching patterns, where the branches are narrow and are twisted with the same handedness (Fig. 1). All colonies of this bacteria grown in similar conditions show the same handedness. Side branches are usually emitted to the convex side of the arced branches (see Refs.  for morphologies and studies of $`𝒞`$ bacteria).

2D chiral branching patterns are also observed in non-living systems. Shapes resembling Sea Horse, or S, are formed during deposition of thin films of fullerene-tetracyanoquinodimethane (C60-TCNQ) or pure TCNQ . Modeling of the system indicates that the apparent curvature of the branches is actually a strong bias in selection of splitting branches, and the branches themselves are not curved. The patterns of bacterial colonies share more resemblance with patterns formed during compression of monolayers of various chiral molecules at the air-water interface , and electro-chemical deposition under a uniform magnetic field . Various modeling attempts indicates that the processes forming these patterns are related to processes of solidification (see Refs.  for the differences between the processes of solidification vs. colonial growth). A different model is required for chiral branching growth of bacterial colonies.

The Communicating Spinors model of Ben-Jacob et al.  is an atomistic (discrete entities) model which describes chiral branching growth of bacterial colonies. The model presented here is not a mean field model of the spinors model. The two models complement each other and highlight different biological features. Each spinor in the spinors model represents a large groups of about $`10^2`$-$`10^4`$ bacteria. This coarse graining makes simulations computationally feasible, but prevents modeling of processes which cannot be averaged over large groups (specifically, single-cell events). Continuous model is more appropriate for this purpose (when interpreting densities as probability densities) . The model presented here can be used to study the effect of mutations and transitions on colonial morphologies, as will be shown elsewhere.

Detailed description of the materials and methods used in the experiments can be found in Ref. . The bacterial colonies are grown on top of agar (semi-solid jelly) with peptone as a nutrient. During colonial growth bacteria are confined to within a layer of lubricating fluid, which is extracted from the agar by the bacteria themselves. Inside the layer the bacteria propel themselves in a random-walk like motion. Optical microscope observations indicate that the length of $`𝒞`$ bacteria is up to 50 times their width. In narrow branches, bacteria are aligned with their neighbors, and their movement is quasi 1D random walk along their long axis.

In order to capture the details of bacterial dynamics we define the bacterial density per angle $`b(\stackrel{}{x},\theta ,t)`$, where $`\stackrel{}{x}𝖱\text{ }^2`$ is a position and $`\theta [0,\pi ]`$ is an angle of orientation. We also use the bacterial density $`B(\stackrel{}{x},t)`$, which is defined as the mean of $`b`$ over $`\theta `$. The equation for the bacterial dynamics is

$$\dot{b}=\left[D_b(B)𝐃_\mathrm{𝟎}(\theta )b\right]_\theta J_\theta +G(b,n)$$
(1)
where $``$ operates in the spatial dimensions, $`D_b`$ is a non-constant diffusion coefficient, $`J_\theta `$ is the flux in the angle dimension which represents changes in bacterial orientation, $`G`$ is a reaction term representing growth and death, and $`n`$ is the nutrient concentration. We take the diffusion matrix $`𝐃_\mathrm{𝟎}\left(\theta \right)=𝐑(\theta )^T\left(\begin{array}{cc}D_d& 0\\ 0& D_l\end{array}\right)𝐑(\theta )`$, where $`𝐑(\theta )`$ is a rotation matrix. $`D_d`$ is a constant coefficient for diffusion in the bacterial direction $`\theta `$ (due to self-propulsion forward and backward) and $`D_lD_d`$ is a constant coefficient for diffusion in the lateral direction, due to fluctuations. Eq. (1) is invariant under rotation in space and translation in $`\theta `$ by the same angle.

It was shown  that bacterial movement in a self-produced layer of fluid can be approximated by a non-linear diffusion, where the diffusion coefficient is proportional to the bacterial density to a power greater than one. The proportion constant is related to the agar dryness through the rate of absorption of the fluid into the agar. Hence we take $`D_b(B)=B^k`$ ($`k2`$), where the proportion constant is included in $`D_d`$ and $`D_l`$. We will assume that bacterial friction with the agar is proportional to the agar dryness, and as the bacterial velocity is inversely proportional to the friction, $`D_d`$ and $`D_l`$ are proportional to the agar dryness to the power $`2`$.

Following  we take a simple form for the growth function, $`G(b,n)=b(n\mu )`$, where $`\mu 0`$ is a rate of conversion into immobile sporulating cells. A more accurate form should have saturation for high values of $`n/B`$, but the linear form is a reasonable approximation for low levels of nutrients, as is the case in the bacterial colonies.

Bacteria change their orientation ($`J_\theta `$) in response to their neighbors orientation. It is convenient to use an auxiliary complex orientation field $`p`$, which is defined as $`p(\stackrel{}{x},t)=\frac{1}{\pi }_0^\pi b(\stackrel{}{x},\theta ,t)e^{i2\theta }𝑑\theta `$. The mean orientation of the bacteria is along the vector $`(\pm \text{Re}\sqrt{p},\text{Im}\sqrt{p})`$.

We found that it is sufficient to track the dynamic of $`B`$ and $`p`$, instead of the full dynamic of $`b`$. Using Fourier expansion of Eq. (1) in $`\theta `$ and taking the first two resulting equations we derive the Orientation Field (OF) model, which includes equations for $`B`$, $`p`$, the nutrient concentration $`n`$ and the density of immobile (sporulating) bacteria $`s`$ (at this stage $`s`$ is used only to record the history of the colonial development):

$`\dot{B}`$ $`=`$ $`\left\{B^k\left[D_1B+2\text{Re}\left(𝐃_\mathrm{𝟐}p\right)\right]\right\}+B(n\mu )`$ (2)
$`\dot{p}`$ $`=`$ $`\left\{B^k\left[𝐃_\mathrm{𝟐}^{}B+D_1p\right]\right\}+p(n\mu )`$ (4)
$`+a(B,|p|)p+\gamma (B^{k2}B,p)ip`$
$`\dot{n}`$ $`=`$ $`D_n^2nnB`$ (5)
$`\dot{s}`$ $`=`$ $`\mu B`$ (6)
where $`n`$ is scaled to units of bacterial mass and $`D_n`$ is its diffusion coefficient. $`a`$ and $`\gamma `$ are real valued functions, a decomposition into orthogonal elements of the derivation of $`_\theta J_\theta `$ (we give here the functions $`a`$ and $`\gamma `$, not the original $`J_\theta `$). Here $`D_1=\frac{D_d+D_l}{2}`$ and $`𝐃_\mathrm{𝟐}=\frac{D_dD_l}{4}𝐌`$ , where $`𝐌=\left(\begin{array}{cc}1& i\\ i& 1\end{array}\right)`$ (from here on we take all constants to be real, unless otherwise stated). Eqs. (2-6) are invariant under a rotation by $`\varphi `$ and a multiplication of $`p`$ by $`e^{i2\varphi }`$.

For the co-alignment function $`a`$ we take

$$a(B,|p|)=4D_\theta +\nu (B|p|)/(B+|p|)$$
(7)
The first term in the RHS results from linear diffusion of $`b`$ in the $`\theta `$ dimension (with $`D_\theta `$ the diffusion coefficient) and the second term is an alignment of non-align bacteria with the mean orientation (with $`\nu `$ being the rate of this co-alignment). The exact form of $`a`$ is not important, as long as $`a\left(|p|\right)`$ has at most one positive root with negative derivative in the range $`[0,B]`$ and $`a\left(|p|\right)`$ in non-positive outside the range $`(0,B)`$.

The rotation function $`\gamma `$ results from the only process that brakes left-right symmetry, the bacterial bias in their tumbling when they are placed in a thin layer of liquid. The bacteria are restrained by their neighbors, and the restraints on rotation are weakest near the boundary of the branches. Due to the non-linear diffusion , the branch’s boundary is identifiable by the divergence of $`B`$ (for the initial conditions specified below). The term $`B^{k1}B`$ is always finite and can be used to identify the boundary (as $`p`$ is bound by $`B`$, it can replace $`B`$ in the coefficient). $`\gamma `$ can be written as a function of $`\mathrm{\Delta }`$, the angle between the branch’s boundary and the bacterial mean orientation:

$`\gamma (\stackrel{}{\alpha },\beta )`$ $`=`$ $`\text{Re}\left[\gamma ^0\left|\beta \stackrel{}{\alpha }\right|\gamma ^0\beta \left(\stackrel{}{\alpha }^T𝐌\stackrel{}{\alpha }\right)/\left|\stackrel{}{\alpha }\right|\right]`$ (8)
$`=`$ $`\left|\beta \stackrel{}{\alpha }\right|\left[2\text{Re}\left(\gamma ^0\right)\mathrm{sin}^2\left(\mathrm{\Delta }\right)+\text{Im}\left(\gamma ^0\right)\mathrm{sin}\left(2\mathrm{\Delta }\right)\right]`$ (9)
where $`\gamma ^0`$ is a complex constant. $`\text{Re}\left(\gamma ^0\right)`$ is a measure of the anti-clockwise bias of the bacteria at the tip of the branch and $`\text{Im}\left(\gamma ^0\right)`$ is a measure of the torque aligning the bacteria in parallel to the boundary. The rotation at the tip of the branch is also restrained by friction with the agar, and from the same reasoning that related $`D_d`$ and $`D_l`$ to the agar dryness, we deduce that $`\text{Re}\left(\gamma ^0\right)`$ is inversely proportional to the agar dryness.

As initial conditions (in a circular 2D geometry), we set $`n`$ to have uniform distribution of level $`n_0`$, $`B`$ to have compact support at the center where it is positive, and the other fields to be zero everywhere. We solve the model numerically using a 2nd order explicit scheme. In order to reduce the implicit anisotropy of the scheme, we use tridiagonal lattice and multiply the bacterial diffusion coefficients by a quenched noise with mean 1. The quenched noise represent inhomogeneities of the agar surface. We show in Fig. 2 that the model can indeed reproduce the microscopic bacterial dynamics and the chiral branching patterns with the local twist.

The Non-Linear Diffusion model (NLD) for colonial growth can reproduce tip-splitting branching patterns  of a related morphotype, $`𝒯`$ morphotype, with shorter bacterial cells. The OF model reduces to the NLD model if the self-propulsion is not primarily along the bacterial long axis ($`D_dD_l`$), if the bacteria do not tend to co-align ($`\nu 4D_\theta `$ or $`\left(\nu 4D_\theta \right)\left(\nu +4D_\theta \right)`$ ), or if bacteria at the tip of the branch tend to rotate too freely ($`\left|\gamma ^0\right|\left(\nu 4D_\theta \right)>\left|𝐃_\mathrm{𝟐}\right|`$). In all these cases $`\left|\text{Re}\left(𝐃_\mathrm{𝟐}p\right)\right|\left|D_1B\right|`$ everywhere and the OF model produce tip-splitting branching patterns.

The response of the simulated growth to initial food concentration and agar dryness is shown in Fig. 3. In agreement with experimental observations , the main effect of the parameters is on the global density of branches, while the curvature of the branches is only weakly affected.

It was shown  that many features of the colonial patterns are explained when food chemotaxis and chemotactic signaling (a chemotactic response to a material emitted by the bacteria themselves) are modeled.

In order to include chemotaxis in our model we bias the spatial flux of $`b`$ in Eq. (1) by adding a chemotaxis term: $`B^kb\zeta _n(n)𝐃_\mathrm{𝟎}(\theta )n`$ (here we take food chemotaxis as an example). $`𝐃_\mathrm{𝟎}(\theta )n`$ is the spatial derivative of the food concentration along the direction of movement of the bacteria. $`\zeta _n(n)𝐃_\mathrm{𝟎}(\theta )n`$ is the derivative as sensed by the bacteria, where $`\zeta _n(n)`$ can be, for example, the “receptor law”  or a constant. The chemotaxis is attractive for positive values of $`\zeta _n(n)`$ and repulsive otherwise.

The bacterial flux due to food chemotaxis translates in the OF model into additional flux terms in Eqs. (2-4). To the flux of $`B`$ in Eq. (2) we add the term $`B^k\zeta _n(n)\left[BD_1𝐈+2\text{Re}\left(p𝐃_\mathrm{𝟐}\right)\right]n`$ (where $`𝐈`$ is the unit matrix) and to the flux of $`p`$ in Eq. (4) we add the term $`B^k\zeta _n(n)\left[𝐃_\mathrm{𝟐}^{}B+pD_1𝐈\right]n`$ .

If a repulsive signaling material $`r`$ is emitted by sporulating cells  then Eqs. (2-4) are affected in a similar manner, with $`\zeta _r(r)`$ replacing $`\zeta _n(n)`$ and $`r`$ replacing $`n`$. $`\zeta _r(r)`$ is negative for a repulsive signal. The equation for the dynamics of the signaling material is

$$\dot{r}=D_r^2r+\mathrm{\Gamma }_rs\lambda _BBr\lambda _rr$$
(10)
where $`\mathrm{\Gamma }_r`$, $`\lambda _B`$, $`\lambda _r`$ are non-negative constants. $`\mathrm{\Gamma }_r`$ is the rate of chemical production by sporulating bacteria, $`\lambda _B`$ is the rate of chemical digestion by bacteria and $`\lambda _r`$ is the rate by which the chemical decompose.

Both food chemotaxis and repulsive chemotactic signaling increase the expansion rate of the colonies. Food chemotaxis does not affect the colonial patterns significantly, while chemotactic signaling does (Fig. 4). For short bacteria, the colonial pattern becomes less ramified, with radial branches and circular global envelope. For long bacteria, the global envelope becomes circular and the branches acquire an outward bias, changing their appearance from an arc-like to a hook-like appearance.

In Fig. 4 we also show the phenomena of global twist . For parameters representing bacteria with intermediate length, repulsive chemotactic signaling can impose global twist on an otherwise tip-splitting pattern. The twist of the branches is relative to the center of the colony, not to the local orientation of the branch. The global nature of the twist is evident by using the “de-chiraling” method presented by Ben-Jacob et al. in Ref.  (see Fig. 4d). In Ref.  Ben-Jacob et al. were able to obtain patterns with a global twist, using the NLD model and applying a rotation operator on the chemotaxis term. As they argue themselves, such operator is inconsistent with the known biological facts and a more detailed model – like the OF model – is required in order to model any type of chirality in the bacterial colonies.

In this letter we present a reaction-diffusion model which accounts for the various morphologies presents by colonies of P. dendritiformis. We focus on chiral features of the colonial patterns, but we are also able to model tip-splitting patterns. We successfully simulated intermediate growth patterns (not shown here). The aim of developing such model is to aid in the study of transition between the two types of growth. Our model is able to handle events of minute density such as a mutation in a single bacterium, and it will be used in following studies of transitions between morphologies and morphotypes.

We thank Inna Brainis for her technical assistance. Presented studies are supported in part by a grant from the Israeli Academy of Sciences grant no. 593/95. IC thanks The Colton Scholarships for their support.


text:
# Effects of dust extinction on optical spectroscopic properties for starburst galaxies in distant clusters

## 1 Introduction

Since Butcher & Oemler (1978) discovered a large fraction of blue galaxies in distant clusters of galaxies at the redshift $`z0.2`$, the origin of the blue color and physical processes closely associated with the formation of the distant blue populations have been extensively discussed by many authors (Dressler & Gunn 1983, 1992; Lavery & Henry 1986; Couch et al. 1994, 1998; Abraham et al. 1996; Barger et al. 1996; Fisher et al. 1998; Morris et al. 1998; Balogh et al. 1999). In particular, Dressler & Gunn (1983, 1992) discovered galaxies with no detectable emission lines and very strong Balmer absorption ones (the so-called E+A galaxies) and concluded that some fraction of distant blue populations are poststarburst galaxies that truncated abruptly their active star formation. One of longstanding and remarkable problems concerning the evolution of galaxies in cluster environment is to understand clearly a physical relationship between these blue populations, whose number fractions rapidly increase with redshift, and passive populations such as elliptical galaxies and S0s. Several attempts have been made to clarify an evolutionary link among star-forming, poststarburst, and passively evolving galaxies observed in distant clusters of galaxies (Couch & Sharples 1987; Abraham et al. 1996; Barger et al. 1996; Morris et al. 1998; Balogh et al. 1999). It is however not so clear whether there is really an evolutionary link among these populations with various spectroscopic properties and what physical process can drive the evolution.

Recent observational studies have found a significant population of possible dusty starburst galaxies in distant clusters (Poggianti et al. 1999: Owen et al. 1999; Smail et al. 1999), which provide a new clue to the evolutionary link among various spectral classes. Poggianti et al. (1999) suggested that galaxies with strong Balmer absorption and relatively modest \[OII\] emission, which are classified as “e(a)” galaxies by Dressler et al. (1999), are dusty starburst galaxies and furthermore discussed that distant clusters have a significant fraction of these dusty starburst galaxies. Smail et al. (1999) found that spectral properties of 5 out of 10 galaxies detected by 1.4 GHz VLA radio observation are classified as poststarburst (a+k/k+a type in Dressler et al. 1999) and considered that the star formation is hidden by dust in them. This finding suggested that dust effects are remarkable even for galaxies that were previously identified as poststarburst ones by optical spectral properties. Although a growing number of observational results on e(a) galaxies with possible dusty starburst have been accumulated, there are only a few theoretical studies addressing the formation and evolution of these e(a) populations in distant clusters.

The purpose of this Letter is to investigate the origin and the nature of e(a) population observed in distant clusters of galaxies by using a one-zone chemical and spectrophotometric evolution model. We particularly investigate whether dusty starburst galaxies can show both strong H$`\delta `$ absorption and relatively modest \[OII\] emission observed in distant e(a) galaxies. We demonstrate that if a younger stellar population formed during a starburst is preferentially obscured by dust in a galaxy, spectral properties of the galaxy during starburst become very similar to those characteristics of e(a) populations. We refer to this behavior of dust extinction as a selective dust extinction, which is originally proposed by Poggianti and Wu (2000) and means that the effect of dust extinction is maximum for the youngest stellar populations among stellar populations with different ages. We furthermore suggest that this selective extinction can be achieved in interacting and merging galaxies in distant clusters. In the followings, the cosmological parameters $`H_0`$ and $`q_0`$ are set to be $`65\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$ and 0.05 respectively, which means that the corresponding present age of the universe is 13.8 Gyr.

## 2 Model

We adopt a one-zone chemical and spectrophotometric evolution model of a disk galaxy with a starburst, and thereby investigate when and how a disk galaxy shows strong H$`\delta `$ absorption line and relatively modest \[OII\] emission one characteristics of e(a) spectra observed in distant clusters of galaxies. Since more details of the adopted model are given in Shioya & Bekki (1998, 2000 in preparation), we only briefly describe the model in the present study. We follow the chemical evolution of galaxies by using the model described in Matteucci & Tornambè (1987) which includes metal-enrichment processes of type Ia and II supernovae (SNIa and SNII). We adopt the Salpeter initial mass function (IMF), $`\varphi (m)m^{1.35}`$, with upper mass limit $`M_{\mathrm{up}}=120M_{}`$ and lower mass limit $`M_{\mathrm{low}}=0.1M_{}`$. We calculate photometric properties of galaxies as follows. The monochromatic flux of a galaxy with age $`T`$, $`F_\lambda (T)`$, is described as

$$F_\lambda (T)=_0^TF_{\mathrm{SSP},\lambda }(Z,Tt)\psi (t)𝑑t,$$
(1)
where $`F_{\mathrm{SSP},\lambda }(Z,Tt)`$ is a monochromatic flux of a single stellar population with age $`Tt`$ and metallicity $`Z`$, and $`\psi (t)`$ is time-dependent star formation rate described later. In the present study, we use the spectral library GISSEL96 which is the latest version of Bruzual & Charlot (1993).

The star formation history of a disk galaxy is characterized by three epochs. The first is the epoch of galaxy formation ($`z_{\mathrm{form}}`$ in redshift and the age of galaxy is 0 Gyr) at which a disk galaxy forms and begins to consume initial interstellar gas by star formation with the moderate rate. In the present study, $`z_{\mathrm{form}}`$ is fixed at 4.5 and therefore the age of galaxies at $`z=0`$ is 12 Gyr. The second is $`z_{\mathrm{sb}}`$ at which starburst begins in the disk. In the present study, $`z_{\mathrm{sb}}`$ is set at 0.4 and the age of galaxies at $`z_{\mathrm{sb}}`$, $`T_{\mathrm{sb}}`$, is 7.64 Gyr. The third is $`z_{\mathrm{end}}`$ at which star formation ceases and is defined as the epoch at which stellar mass fraction becomes 0.95 in our models. The age of galaxies at $`z_{\mathrm{end}}`$, $`T_{\mathrm{end}}`$, is 8.56 Gyr. In the following, we mainly use the age of galaxies to describe their evolution. Throughout the evolution of disk galaxies, the star formation rate is assumed to be proportional to gas mass fraction ($`f_g`$) of galaxies;

$$\psi (t)=kf_g,$$
(2)
where $`k`$ is a parameter which controls the star formation rate. This parameter $`k`$ is given as follows:

$$k=\{\begin{array}{ccc}k_{\mathrm{disk}}\hfill & & \mathrm{for}\mathrm{\hspace{0.33em}\hspace{0.33em}0}T<T_{\mathrm{sb}},\hfill \\ k_{\mathrm{sb}}\hfill & & \mathrm{for}T_{\mathrm{sb}}T<T_{\mathrm{end}},\hfill \\ 0\hfill & & \mathrm{for}T_{\mathrm{end}}T.\hfill \end{array}$$
(3)
In the following, we refer the first, second, and third phases as the prestarburst, the starburst, and the poststarburst, respectively. In the present study, we set the value of $`k_{\mathrm{disk}}`$ as 0.225 $`\mathrm{Gyr}^1`$, which corresponds to a plausible star formation rate for Sb disks (e.g., Arimoto, Yoshii, & Takahara 1992). $`k_{\mathrm{sb}}`$ is considered to be 10 times larger than $`k_{\mathrm{disk}}`$ ($`2.25\mathrm{Gyr}^1`$), which is consistent with observational results on starburst galaxies (e.g., Planesas et al. 1997). Figure 1 shows the star formation history of our model.

We try to understand the role of the selective dust extinction in the formation of e(a) galaxies by comparing the following two models with each other; a model with no dust extinction (from now on referred to as the $`ND`$ model) and the selective dust extinction model (the $`SD`$ model). In the $`ND`$ model, dust effects on photometric and spectroscopic properties are completely neglected; i.e., $`A_V=0`$ at all times. In the $`SD`$ model, we assume that during starburst the value of $`A_V`$ depends on the age of stellar population; i.e.,

$$A_V=A_{V,0}\mathrm{exp}\{(Tt)/\tau \}\mathrm{for}T_{\mathrm{sb}}<T<T_{\mathrm{end}},$$
(4)
where $`A_{V,0}`$ and $`\tau `$ are parameters controlling the degree of extinction. We here set $`A_{V,0}=5`$ mag and $`\tau =1.0\times 10^6`$ yr.

Based on the monochromatic flux derived in the equation (1) and the value of $`A_V`$ from the equation (4) for each age of a galaxy $`T`$, we calculated the SED of a galaxy corrected by dust extinction using the extinction law derived by Cardelli et al. (1989) and adopting the so-called screen model. For deriving the fluxes for various gaseous emission lines (H$`\delta `$ and \[OII\]) in dusty starburst galaxies, we first calculate the number of Lyman continuum photons, $`N_{\mathrm{Ly}}`$, by using the SED that are not modified by dust extinction. If all the Lyman continuum photons are used for ionizing the surrounding gas, the luminosity of H$`\beta `$ is calculated according to the following formula;

$$L(\mathrm{H}\beta )(\mathrm{erg}\mathrm{s}^1)=4.76\times 10^{13}N_{\mathrm{Ly}}(\mathrm{s}^1)$$
(5)
(Leitherer & Heckman 1995). To calculate luminosities of other emission lines, e.g., \[OII\] and H$`\delta `$, we use the relative luminosity to H$`\beta `$ luminosity tabulated in PEGASE (Fioc & Rocca-Volmerange 1997) which is calculated for the set of electron temperature of 10000 K and electron density of 1 cm<sup>-3</sup>. Thus the SED derived in the present study consists of stellar continuum and gaseous emission.

## 3 Results

Figure 2 shows the time evolution of EW(H$`\delta `$) and EW(\[OII\]) for the $`ND`$ model and the $`SD`$ one. First, we describe the evolution for the $`ND`$ model. The equivalent width of \[OII\] decreases gradually before the starburst begins (i.e., at $`T<7.64\mathrm{Gyr}`$). At $`T=T_{\mathrm{sb}}`$, the equivalent width of \[OII\] increases rapidly following the rapid increase of star formation rate due to starburst and then decreases gradually with the decreasing star formation rate. During the starburst, the equivalent width of \[OII\] is larger than 40 Å. After star formation ceases at $`T=T_{\mathrm{end}}`$ and massive stars formed during starburst die out, the equivalent width of \[OII\] rapidly decreases and becomes smaller than 5 Å. The strength of equivalent width of H$`\delta `$ evolves as follows. During the prestarburst evolution phase, the H$`\delta `$ line is observed as an absorption one and its equivalent width reaches the local maxima at around 3 Gyr ($`z1.5`$). The strength of EW(H$`\delta `$) gradually decreases, and becomes about 3 Å just before starburst begins ($`T=7.64\mathrm{Gyr}`$, $`z=0.4`$). When starburst begins, the H$`\delta `$ line changes from an absorption line to an emission one and its equivalent width reaches $`5`$ Å. $`\mathrm{A}\mathrm{few}\times 10^8\mathrm{yr}`$ after the starburst begins, H$`\delta `$ line changes into absorption line again. Since the flux of H$`\delta `$ emission line becomes negligibly small, the equivalent width of H$`\delta `$ reaches 8 Å when starburst ceases. As a result of passive evolution of galaxies after starburst, the equivalent width of H$`\delta `$ decreases gradually: it becomes smaller than 7 Å at 8.8 Gyr and 3 Å at 11 Gyr.

Next we describe the evolution of the $`SD`$ model. Since we assume that the selective extinction affects spectroscopic properties of galaxies only in starburst phase, there is no difference in the evolution of EW(\[OII\]) and EW(H$`\delta `$) between the $`SD`$ model and the $`ND`$ one in the prestarburst phase and the poststarburst one. We therefore pay our attention to the evolution of EW(\[OII\]) and EW(H$`\delta `$) during starburst. In the $`SD`$ model, the equivalent width of \[OII\] does not increase (but increases in the $`ND`$ model) and is below 40 Å during starburst. The equivalent width of H$`\delta `$ becomes rapidly small just after starburst begins, although the H$`\delta `$ line can be still observed as an absorption line. A few $`\times 10^8`$ years after the starburst begins, the equivalent width of H$`\delta `$ becomes larger than 4 Å. The reason why the equivalent width of H$`\delta `$ of the $`SD`$ model is always larger than that of the $`ND`$ model is that the flux of H$`\delta `$ emission line of the $`SD`$ model is smaller than that of $`ND`$ model.

Figure 3 shows the evolution of galaxies on EW(\[OII\])-EW(H$`\delta `$) plane for each model. The criteria of spectroscopic classification determined by Dressler et al. (1999) are also superimposed on it. The $`ND`$ model is located within the e(b) region during starburst and is very rapidly shifted to the a+k one just when starburst is ended. Although, the locus of the $`ND`$ model passes through the e(a) region, the crossing time is $`4\times 10^6`$ years which is about 0.5 % of the duration of the starburst. On the other hand, the $`SD`$ model can successfully show the spectroscopic properties of e(a) galaxies, namely, strong H$`\delta `$ absorption line \[EW (H$`\delta `$) $`>4`$ Å\] and relatively modest \[OII\] emission line \[EW(\[OII\]) $`>`$ 5 Å\], during the most of the starburst phase (0.7 Gyr). This result accordingly suggests that a selective extinction plays a vital role in the formation of galaxies with e(a) spectra. These results thus confirm the early suggestion by Poggianti & Wu (2000) that if the youngest stellar population is the most heavily obscured by dust in a dusty starburst galaxy, the e(a) spectra can be achieved in the galaxy.

As is described above, we have confirmed that the $`SD`$ model shows the weaker \[OII\] emission and stronger H$`\delta `$ absorption compared with the $`ND`$ model. The physical reason for the successful reproduction of the $`SD`$ model can be understood in terms of the difference in the influence of dust extinction between continuum and emission lines. The emission lines come from the ionized gas around youngest stellar populations which are most heavily obscured by dust. On the other hand, most of the continuum come from the stellar populations that are older and less obscured than younger ones dominating the ionizing photons. The flux of emission lines consequently is more greatly affected by dust extinction than that of stellar continuum. Thus the $`SD`$ model can show the smaller equivalent width of emission lines without changing so greatly the flux of continuum.

## 4 Discussion and conclusions

Poggianti & Wu (2000) first discussed that if the dust extinction for a stellar population decreases with the stellar age, e(a) spectra can be obtained for a starburst galaxy. They furthermore suggested that if the location and thickness of dust patches depend on the age of the embedded stellar populations in a starburst galaxy, the effects of the above selective dust extinction become remarkable for the galaxy. Then, when and how the proposed difference in the degree of dust extinction between stellar populations with different ages is possible during galactic evolution in distant clusters? We here suggest that the difference of dust extinction between the central region and the outer one in a galaxy is a main cause for the above age-dependent extinction. To be more specific, since younger stellar populations formed by secondary starburst in the central region with a larger amount of dusty interstellar gas are more heavily obscured by dust than the outer old components, the age-dependent extinction can be achieved. This radial dependence of dust extinction is considerably reasonable and realistic, considering that secondary dusty starburst occurs preferentially in the central part of a galaxy owing to efficient inward gas transfer driven by non-axisymmetric structure (such as stellar bars) and galaxy interaction and merging. Numerical simulations on gaseous and stellar distribution in merging disk galaxies with dusty starburst furthermore have demonstrated that the central young stellar component formed by nuclear starburst is more heavily obscured by dusty gas than the outer old component initially located in merger progenitor disks (Bekki, Shioya, & Tanaka 1999). Thus we strongly suggest that galaxies with the central dusty starburst are more likely to show the selective dust extinction. Future observational studies on the radial dependence of photometric and spectroscopic properties for e(a) galaxies, which can reveal the detailed distribution of young dusty population and that of old one with less extinction, will assess the validity of the above idea for the origin of the selective dust extinction (i.e., the age-dependent extinction).

Then what physical process is closely associated with the radial dependence of dust extinction suggested above and thus with the formation of e(a) galaxies in distant clusters ? We here suggest that minor and unequal-mass galaxy merging, which is different from major merging in that it can leave disk systems even after merging, is one of possible candidates that can reproduce fundamental properties of e(a) galaxies. Mihos & Hernquist (1994) has demonstrated that minor merging can excite non-axisymmetric structure in a late-type disk, trigger nuclear starburst, and transform the disk into an early-type disk such as Sa and S0. Unequal-mass merging between two disks is furthermore demonstrated to create a disk with the structure and kinematics strikingly similar to those of S0s (Barnes 1994), and furthermore Bekki (1998) found that such merging can trigger nuclear starburst, form a very gas-poor disk, and create a remarkable bulge characteristic of typical S0s. These merging therefore not only contribute to the formation of dusty nuclear starburst and thus to the radial dependence of dust extinction but also provide an evolutionary link between late-type disks, e(a) galaxies, and early-type disks (S0) with the number fraction observed to decrease with redshift (e.g., Dressler et al. 1996). An early-type disk galaxy formed by minor and unequal-mass merging experiences the post-starburst phase with k+a/a+k spectra after starburst and shows a remarkable bulge with young stellar populations and thus with bluer colors (Bekki 1998). Accordingly one of observational tests for assessing the relative importance of minor and unequal-mass merging in the fraction of e(a) galaxies with the selective dust extinction is to investigate the number fraction of disk galaxies with k+a/a+k spectra and prominent blue bulges in distant clusters.

We are grateful to the anonymous referee for valuable comments, which contribute to improve the present paper. Y.S. thanks to the Japan Society for Promotion of Science (JSPS) Research Fellowships for Young Scientist.


text:
# Electron transport in nanotube–molecular wire hybrids

## Abstract

We study contact effects on electron transport across a molecular wire sandwiched between two semi-infinite (carbon) nanotube leads as a model for nanoelectrodes. Employing the Landauer scattering matrix approach we find that the conductance is very sensitive to parameters such as the coupling strength and geometry of the contact. The conductance exhibits markedly different behavior in the two limiting scenarios of single contact and multiple contacts between the molecular wire and the nanotube interfacial atoms. In contrast to a single contact the multiple-contact configuration acts as a filter selecting single transport channels. It exhibits a scaling law for the conductance as a function of coupling strength and tube diameter. We also observe an unusual narrow-to-broad-to-narrow behavior of conductance resonances upon decreasing the coupling.

PACS numbers: 73.50.-h, 73.61.Wp 85.65.+h

The ‘top-down’ miniaturization of electronic devices to micro- or nanometer-size has triggered formidable progress in the research of mesoscopic systems for more than a decade now. Such endeavor has given rise to a variety of nontrivial quantum effects in solid state physics at these length scales providing an active field of interplay between fundamental and applied research.

An arena closely related to mesoscopics is molecular electronics which combines molecular and solid state physics. Driven by possible technological applications and owing to recent experimental breakthrough, this area has evolved to a rapidly developing field on its own posing new challenges to theory. Molecular electronics is based on the ‘bottom-up’ manufacture philosophy whose underlying principle is to use molecules or supramolecular structures as (reproducible) circuit elements. Hence, employing the diverse electronic properties of molecular complexes and the capacity of synthetic (bio)chemistry, it may be able to circumvent inherent limitations and imperfections of conventional semiconductor device fabrication techniques at nanometer scales.

Although the original idea is quite old, significant progress has only been demonstrated experimentally in recent years. Owing to the advances in self-assembly techniques , end-group modifications , scanning probe  and break-junction techniques , atomic-scale control and positioning of single molecules and their assemblies become possible. First electron transport measurements through molecular complexes between metallic electrodes have been reported. Proposals and studies of molecular wires range from ‘simple’ molecules  to DNA strands . In a parallel development the use of carbon nanotube networks has been the focus of intense experimental and theoretical activity as another promising direction for building blocks of molecular circuits .

Albeit a molecular device is typically divided into three parts, the donor and acceptor electrodes and the molecular compound serving as a bridge, it is clear that to understand conductance measurements an account of the system as a whole is required. This is intuitively plausible when looking at the problem as an electron transfer process. The coupling of the molecular complex to the environment as well as its intrinsic ability to convey charge are equally important factors. Hence, whereas the molecular character has been the main focus , the precise nature of the contact and its implications has also become a topic of investigation .

Yet the electrodes are usually formed from bulk material. In contrast, we take here the viewpoint that the electrodes in the vicinity of the molecular interface can be mesoscopic themselves. To be more specific, we focus on carbon tubules as suitable candidates for such nanoelectrodes. Carbon nanotubes are known to exhibit a wealth of properties depending on their diameter ($`nm`$), chirality (orientation of graphene roll up), and whether they consist of a single cylindrical surface (single-wall) or more (multi-wall) . On the one hand, first experimental attempts to build nanotube-molecule-nanotube hybrids are on their way . On the other hand, carbon nanotubes are utilized as scanning probe tips to study molecular structures . This represents a related setup where contact effects of a molecule-nanoelectrode junction play a key role.

In the present study, we address the influence of the molecular wire-electrode contact on the conductance for the class of systems where the structure of electrodes plays an important role. For mesoscopic leads with reduced dimensionality lateral to the current direction, it makes sense to discriminate between electron transport channels, e.g. carbon nanotubes support up to two channels for electrons with energy around the equilibrium Fermi energy. Evidently, for such low-dimensional transport the geometry of the contact should crucially determine the measured conductance. We find that electron transport shares distinct properties depending on the number and strength of contacts between the molecular bridge and the interface as well as on the symmetry of the channel wavefunctions transverse to transport. We demonstrate that single contacts give rise to complex conductance spectra exhibiting quantum features of both the molecule and the electrodes; multiple contacts provide a mechanism for transport channel selection, leading to a scaling law for the conductance and allowing for its control. Channel selection also highlights the role of molecular resonant states by suppressing details assigned to the electrodes. Such information may be used as a guideline for systematic chemical synthesis or complementary experimental analysis.

We shall now specify more the system we have in mind. The electrodes are open-ended single-wall carbon nanotubes described by a parameterized tight-binding Hamiltonian with a $`\pi `$-electron per atom which represents a good approximation to the carbon nanotubes electronic band structure , whereas the molecular system sandwiched in between is modeled by a homogeneous tight-binding chain. The latter choice reflects our aim to provide a qualitative understanding of the phenomenology of the physical problem as opposed to more quantitative quantum-chemical methods that take into account the precise structure and properties of the molecular bridge.

The electronic Hamiltonian of the full system, including the left (L) and right (R) tube (see Fig. 1), reads

$`H`$ $`=`$ $`H_{\mathrm{tubes}}+H_{\mathrm{wire}}+H_{\mathrm{coupling}}`$ (1)
$`=`$ $`{\displaystyle \underset{\alpha =\mathrm{L},\mathrm{R},\mathrm{wire}}{}}{\displaystyle \underset{n_\alpha ^{},n_\alpha ^{}}{}}\left({\displaystyle \frac{\epsilon _{n_\alpha }^\alpha }{2}}\delta _{n_\alpha ^{},n_\alpha ^{}}+\gamma _{n_\alpha ^{},n_\alpha ^{}}^\alpha \right)|n_\alpha ^{}n_\alpha ^{}|`$ (2)
$`+`$ $`{\displaystyle \underset{m_\mathrm{L}}{}}\mathrm{\Gamma }|m_\mathrm{L}n_{\mathrm{wire}}=1\left|+{\displaystyle \underset{m_\mathrm{R}}{}}\mathrm{\Gamma }\right|m_\mathrm{R}n_{\mathrm{wire}}=N|+\mathrm{H}.\mathrm{c}..`$ (3)
Here, $`\gamma ^{\mathrm{L},\mathrm{R}}`$($`=2.66`$ eV), $`\gamma ^{\mathrm{wire}}`$, and $`\mathrm{\Gamma }`$ are the hopping matrix elements between atoms of the carbon tubule leads, molecular bridge, and the bridge/lead interface, respectively. They are non-zero only for nearest neighbors. In Eq. (1), $`\epsilon ^{\mathrm{wire}}`$ is the on-site or orbital energy of each of the $`n_{\mathrm{wire}}=1,\mathrm{},N`$ chain-atoms relative to that of the leads, $`\epsilon ^{\mathrm{L},\mathrm{R}}`$, which is fixed to zero. Summations over $`m_\mathrm{L}`$ and $`m_\mathrm{R}`$ run over interfacial end-atoms of the leads. In general, there are $`M`$ such atomic positions, depending on the perimeter of the tubes, and $`1M_\mathrm{c}M`$ hybridization contacts. We also compare with a square lattice model of mesoscopic electrodes with nearest-neighbor interactions ($`\gamma ^{\mathrm{L},\mathrm{R}}=1`$eV) and periodic boundary conditions, which delivers additional insight.

In what follows, we use the Landauer theory  which relates the conductance of a system to an independent-electron scattering problem  and describes unique quantum effects in mesoscopic systems . The electron wavefunction is assumed to extend coherently across the device and the two-terminal, linear-response conductance at zero temperature reads

$$G(E_\mathrm{F})=2(e^2/h)T(E_\mathrm{F}).$$
(4)
The factor two accounts for spin degeneracy, and $`T(E_\mathrm{F})`$ is the total transmittance for injected electrons with Fermi energy $`E_\mathrm{F}`$. The transmission function is given by $`T(E)=_{\nu \nu ^{}}|S_{\nu \nu ^{}}(E)|^2`$, where $`\nu ,\nu ^{}`$ are quantum numbers labeling open channels for transport which belong to mutually exclusive leads, in our case the two semi-infinite perfect nanotubes. The molecular system attached acts as a scatterer, and $`S`$ is the corresponding quantum-mechanical scattering matrix.

For the numerical calculation of the central quantity, $`T(E)`$, we use a general scattering technique which has been recently formulated for studies of the giant magnetoresistance . Application of the method to phonon transport across disordered interfaces  and to electric conductance in multi-wall carbon nanotubes  reveals an efficient algorithm for calculating the Green function for arbitrary tight-binding Hamiltonians and, hence, the $`S`$-matrix . The computational scheme comprises two essential steps: first, the calculation of an effective (renormalized) interaction between the electrodes by projecting out the degrees of freedom of the scatterer, and second, the computation of the unperturbed electrodes Green function. Then, one uses the Dyson equation to express the Green function of the composite system (leads plus scatterer). We followed the suggestion discussed in  for the implementation and computation of the exact Green function of the tubules.

Gross properties of the conductance spectrum of the system can be understood by looking at the two extreme cases of a single interfacial contact, $`M_\mathrm{c}=1`$ (single contact-SC), and multiple contacts, $`M_\mathrm{c}=M`$ (MC). In the SC-scenario all open channels contribute to the transmission, i.e. $`S_{\nu \nu ^{}}(E)`$ is non-zero for any $`\nu ,\nu ^{}`$. For the case of the molecular wire bridging two square-lattice tubes, depicted in the upper left panel of Fig. 2, the conductance bears some of the properties of current flowing through a one-dimensional constriction . In particular, the conductance shows resonances of quantum unit ($`2e^2/h`$) height at eigenenergies of the isolated molecular chain (indicated as triangles). They arise because of back reflections at the molecular interface. However, for carbon tubules leads (upper right panel of Fig. 2) we observe additional structure in the conductance spectra. Preliminary results suggest that distinctive features such as antiresonances are signatures of van Hove singularities in the carbon tubules band structure .

In contrast, the MC-configuration acts as a channel filter resulting in a profoundly different behavior: The conductance vanishes for part of the spectrum as shown in the lower panels of Fig. 2. The complicated conductance spectrum for a SC carbon tubule-molecule configuration turns into a regular sequence of peaks at eigenenergies of the isolated molecular wire (marked as diamonds). Further analysis of the $`S`$-matrix elements revealed that only wavefunctions of the tubes without modulation along the cross-section circumference allow for transport, thereby, yielding zero conductance when such channels are not available. The filtering is a consequence of a sum rule that determines the transmission of each open channel which may be viewed as the overlap $`n_{\mathrm{wire}}|H|n_{\mathrm{L},\mathrm{R}}`$ (see Eq. (1), due to the nature of coupling only the transverse profile is important). The overlap is also related to the spectral density . For the square lattice model tubes e.g. only the channel with zero ‘transverse’ momentum gives a non-zero sum. From symmetry considerations it is clear that the channel selection is a generic feature of cylindrical electrodes. Moreover, we found that channel filtering approximately prevails also for non-cylindrical, mesoscopic electrodes with lateral confinement. More generally, multiple contacts allow for control of low-dimensional transport via channel selection.

An additional particular feature of the MC configuration is that the conductance conforms to a scaling law. Two MC hybrid–structures differing in tube diameters $`D`$ and contact strengths $`\mathrm{\Gamma }`$, but conserving the product $`\mathrm{\Gamma }^2D`$, exhibit the same conductance profile (Fig. 3). This is a mere contact effect related to the symmetry of the contributing channel wavefunction and should hold for any effective coupling with the form considered here. For single-channel transport it follows that the conductance is proportional to $`|n_{\mathrm{wire}}|H|n_{\mathrm{L},\mathrm{R}}|^2`$. Then one readily arrives at the exact form of the scaling law taking into account that the transverse profile of the contributing channel wavefunction has no nodes and is normalized.

The ‘intermediate contact’ case, i.e. $`1<M_c<M`$, exhibits much richer behavior but can be understood with the above arguments. For completeness we discuss a specific example which once more illustrates the importance of the interfacial coupling for molecular systems bridging nanoelectrodes and further supports channel selection. In Fig. 4 typical conductance spectra for zigzag carbon nanotube electrodes are shown. To interpret these results we note that transport usually takes place at $`E1`$eV around the Fermi energy $`E_F=0`$. For this part of the spectrum we notice a complete suppression of conductance for $`M_c=3`$, owing to contact ‘dimensionality’ . The origin of this effect derives from metallic zigzag nanotubes supporting two degenerate transport channels in this energy region with wavefunction symmetries such that the wire/tube overlap gives a zero contribution for $`M_c=3n`$ as depicted and non-zero otherwise .

Finally, we examine the dependence of the conductance on the coupling strength $`\mathrm{\Gamma }`$. It is expected that as the coupling is reduced, resonances are resolved with decreasing width-to-height ratio. This is indeed visible when comparing the middle and lower panels of Fig. 5 for both single (left) and multiple contact (right). However, the upper panels of Fig. 5 reveal a counterintuitive behavior. By decreasing $`\mathrm{\Gamma }`$ from top to bottom in Fig. 5 we observe an intermediate broadening of the resonance widths, accompanied by a clear-cut change in the position and number of resonances picked up. There is a crossover from $`N2N`$ resonances, if all the wire resonant states (marked as diamonds and triangles in Fig. 5) fall into the tube energy band. The position of resonances is suggestive for a possible explanation of this effect. It points to considering the wire end-atoms as part of the leads due to strong binding with the transition driven by the interaction strength. Similar binding effects, which determine the effective length (or the nature) of the molecular bridge, have been discussed in Ref.  and are quite sensitive to multiple contacts as demonstrated.

In additional simulations , we have considered a dimerized chain as a model for a molecular wire as suggested in  and more complex artificial molecular compounds, which nevertheless do not show distinct behavior from what we have discussed.

We have not discussed electron transport for finite temperatures and bias voltages, since our objective has been to point out generic geometry-induced features linked to the contact resistance of molecular wires. In that case the calculation of $`I`$-$`V`$ characteristics requires generalized methods with the Landauer conductance delivering the initial insight. We further point out that in the scattering approach used it is implicitly assumed that no inelastic scattering or other phase-breaking mechanisms occur. For instance, coupling to molecular vibrational modes is neglected and, therefore, studies relying on the independent-electron picture may be only indicative of the underlying conduction mechanisms. We have also assumed that the molecular system is strongly coupled to the electrodes so that charging effects are not important.

In summary, we have studied contact effects on the conductance of a hybrid built by a molecular wire coupled to low-dimensional leads. We have pointed out immediate consequences of the contact geometry and dimensionality on electron transport across such systems. These include channel selection, conductance spectra tuning, and the existence of a scaling law $`G=G(\mathrm{\Gamma }^2M)`$. We also demonstrated that a square lattice tube exhibits most of the above features and can be used as a tractable analytic model. Replacing such a model for electrodes by natural and realistic candidates of molecular electronic circuits, carbon nanotubes, adds to the conductance profile richer structure. The latter needs further investigation to relate it to the tubules band structure and details depending on the tubule chirality. We hope that our observations may drive systematic chemical synthesis based on geometrical factors and stimulate corresponding experimental analysis. Especially for multi-wall nanotube leads, an axial magnetic field of reasonable magnitude can be applied which modulates the symmetry of the tube states and, hence, can act as an external tuning parameter.

We acknowledge fruitful discussions with A. Bachtold, V.I. Fal’ko, and P. Fulde. G.C. was partially supported by the EU under contract TMR-ERBFMRXCT980180.


text:
# The Ly𝛼 Forest Flux Distribution at 𝑧∼5.2 and the Evolution of the Ionizing Background

## 1 INTRODUCTION

The theory of the Ly$`\alpha `$ forest based on structure formation in Cold Dark Matter models has been highly successful in accounting for the observable properties of the Ly$`\alpha `$ forest, in particular the flux distribution and the power spectrum (see Croft et al. 1999, McDonald et al. 2000a, and references therein). As new quasars at ever higher redshift are discovered, the comparison of these observable quantities with the theory should be extended. Because the theory predicts the evolution of the spatial distribution of intergalactic gas, it is now possible to infer the evolution of the ionizing background intensity from the observed evolution of the mean flux decrement. Furthermore, as the epoch of reionization is approached, fluctuations in the background intensity due to the discreteness of the sources should increase, affecting the Ly$`\alpha `$ forest flux distribution.

In this paper, we use the observations of three quasars at $`z5`$ (Songaila et al. 1999, hereafter SHCM; Stern et al. 2000; Fan et al. 2000) to infer the evolution in the hydrogen ionization rate, and to test if the flux distribution is consistent with theoretical predictions. In §2 we present predictions of the flux distribution at high redshift, using numerical simulations whose resolution and box size are shown to be sufficient to have produced convergent results. In §3 we compare the predictions to the observations to measure the evolution of the ionization rate and search for differences in the predicted and observed flux distribution.

## 2 NUMERICAL RESULTS FOR THE FLUX DISTRIBUTION

If we ignore the effects of peculiar velocity of the gas and of thermal broadening, then the optical depth to Ly$`\alpha `$ scattering, $`\tau `$, at any given point of an observed spectrum is determined by the density of gas at the corresponding point in space, $`\rho _{\mathrm{gas}}\mathrm{\Delta }(\mathrm{\Omega }_b\rho _{\mathrm{crit}})`$, and by the gas temperature, which we assume to follow the relation $`T=T_0\mathrm{\Delta }^{\gamma 1}`$ (see Hui & Gnedin 1997):

$$\tau =\tau _0\frac{\left(1+z\right)^6\left(\mathrm{\Omega }_bh^2\right)^2}{T_0^{0.7}H(z)\mathrm{\Gamma }_{12}(z)}\mathrm{\Delta }^\beta .$$
(1)
Here, $`\tau _0`$ is a constant, $`H(z)H_0\mathrm{\Omega }_m^{1/2}(1+z)^{3/2}`$ is the Hubble constant at redshift $`z`$, $`H_0`$ is the present Hubble constant, $`h=H_0/(100\mathrm{km}\mathrm{s}^1/\mathrm{Mpc})`$, $`\mathrm{\Omega }_m`$ is the fraction of the critical density in pressureless matter, $`\mathrm{\Gamma }_{12}(z)`$ is the photoionization rate of hydrogen in units of $`10^{12}\mathrm{s}^1`$, and $`\beta =20.7(\gamma 1)`$ (where we use $`\alpha (T)T^{0.7}`$ for the recombination coefficient). The dependence of the distribution and average of the transmitted flux, $`F=\mathrm{exp}(\tau )`$ (which is the observed quantity), on the cosmological model can be understood from equation (1). Thus, as long as the distribution of $`\mathrm{\Delta }`$ and the $`\mathrm{\Delta }T`$ relation are sufficiently well understood, the observed mean transmitted flux of the Ly$`\alpha `$ forest can be used to determine the parameter combination $`(\mathrm{\Omega }_bh^2)^2/H(z)/\mathrm{\Gamma }_{12}(z)`$ in equation (1).

Knowledge of the distribution of $`\mathrm{\Delta }`$ and $`T`$ needs to come from the theory, with the use of hydrodynamic simulations. This theory can be tested by checking its predictions for the detailed form of the distribution of $`F`$ and its power spectrum against the observations, and this comparison has shown that the Cold Dark Matter model with a cosmological constant that is supported by a large body of observational evidence in galaxy clustering and CMB fluctuations (Primack 2000, and references therein) is in good agreement with the observed Ly$`\alpha `$ forest (Rauch et al. 1997, Croft et al. 1999, McDonald et al. 2000a, and references therein). At the same time, the gas temperature over the relevant range of densities has been measured from the absorption line shapes (Ricotti, Gnedin, & Shull 2000; Schaye et al. 2000; McDonald et al. 2000b), although the systematic differences between these studies are still substantial.

We focus in this paper on the flux distribution at high redshift. This particular property of the Ly$`\alpha `$ forest is affected only by the low-density intergalactic gas, given the large mean flux decrement (the reason is that high density gas causes absorption only in highly saturated regions of the spectrum). We first show that this distribution is predicted without significant uncertainties due to the finite resolution and box size of the simulations. We use the Hydro-PM (HPM) code of Gnedin & Hui (1998), comparing it first to the fully hydrodynamic simulation in Miralda-Escudé et al. (1996) to test that the approximate nature of the Hydro-PM method does not result in any significant differences in the predicted flux distribution at high redshift. When we create simulated spectra, we assume a $`\mathrm{\Delta }T`$ relation given by $`T_0=2\times 10^4`$ K and $`\gamma 1=0`$, which is consistent with the observational determinations Ricotti et al. (1999) and McDonald et al. (2000b) at $`z4`$, although Schaye et al. (2000) find a lower temperature (for the calculation of the pressure term in the HPM simulations we assume $`T=0`$ before reionization at $`z=7`$, and then interpolate between $`T_0=25000`$, $`\gamma 1=0.0`$ at $`z=7`$ and $`T_0=19000`$, $`\gamma 1=0.2`$ at $`z=3.9`$; however, we find that the results are almost completely insensitive to these assumptions, probably because the pressure does not have much influence in the voids). The cosmological model is $`\mathrm{\Lambda }`$CDM, with $`\mathrm{\Omega }_m=0.4`$, $`\mathrm{\Omega }_bh^2=0.015`$ (note, however, that the results we give for $`\mathrm{\Gamma }_{12}`$ are rescaled to $`\mathrm{\Omega }_bh^2=0.02`$), $`\mathrm{\Omega }_\mathrm{\Lambda }=0.6`$, $`h=0.65`$, and $`n=0.95`$, where $`\mathrm{\Omega }_\mathrm{\Lambda }`$ is the fraction of the critical density in cosmological constant and $`n`$ is the usual logarithmic slope of the primordial power spectrum. Simulated spectra are constructed in which the transmitted flux is binned into 2 Å pixels, to match the observed spectrum in Fan et al. (2000) (we use this pixel size throughout the paper). In Figure 1, we show the flux distribution predicted for the full-hydro simulation as the heavy solid line. The light (long-dashed, dotted, solid, dashed) lines show HPM simulations with box size $`L=`$(17.8, 8.9, 4.4, 4.4), and a total number of particles $`N=(256^3,128^3,128^3,256^3)`$, respectively. All the results are at $`z=4`$, but using the value of $`\mathrm{\Gamma }_{12}`$ that gives a mean flux decrement from the hydrodynamic simulation of $`\overline{F}=0.1`$, which matches the observations at $`z5`$. We keep $`\mathrm{\Gamma }_{12}`$ (rather than $`\overline{F}`$) constant for all the simulations, so that we can see the difference in the predicted $`\overline{F}`$.

The result for the fully hydrodynamic simulation and for the HPM simulation that has almost the same box size and resolution ($`L=8.9`$ and $`N=256^3`$, not plotted) are essentially identical. This justifies using the HPM simulations to examine the effects of changing the box size and resolution and to compute the flux distribution for comparison with observations at different redshifts. The two simulations with the smallest box size, with resolutions $`L/N^{1/3}=35h^1\mathrm{Kpc}`$ and $`17h^1\mathrm{Kpc}`$ (solid and dashed lines, respectively), show that simulations with inter-particle separation of $`35h^1\mathrm{Kpc}`$ have already achieved the required resolution; however, degrading the resolution to $`L/N^{1/3}=69h^1\mathrm{Kpc}`$ causes a change in $`\overline{F}`$ of $`0.02`$, similar to the observational error bars of $`\overline{F}`$. Comparison of the long-dashed and dotted (or the heavy and light solid) lines shows that the effect of increasing the box size beyond $`5h^1\mathrm{Mpc}`$ is small. We conclude that simulations with boxes $`L=8.9h^1\mathrm{Mpc}`$ and number of particles $`N=256^3`$ are sufficient for our purpose to have achieved convergence. We find similar convergence in the HPM simulations at $`z=5.2`$. We also find that different random realizations of the initial conditions of the simulation of this box size and resolution change $`\overline{F}`$ by $`\pm 0.01`$ (small compared to the observational error).

## 3 RESULTS

We now examine the flux distribution at high redshift predicted from the simulations to do a direct comparison to the observational constraints. We use the HPM simulation with $`L=8.9h^1\mathrm{Mpc}`$ and $`N=256^3`$, and create spectra binned into 2 Å pixels.

The curves in Figure 2(a) show the cumulative flux distribution at $`z=5.2`$, for ten values of $`\mathrm{\Gamma }_{12}`$, corresponding to different values of $`\overline{F}`$, indicated in the Figure. This can be compared to the spectrum in Figure 3 of Fan et al. (2000), from a quasar at $`z=5.8`$ (where the mean redshift of the observed Ly$`\alpha `$ forest, in the wavelength range 7140-7956 Å used by Fan et al. 2000 to compute $`\overline{F}`$, is $`5.2`$). In order to match the observed $`\overline{F}=0.10\pm 0.02`$, we need $`\mathrm{\Gamma }_{12}0.17\pm 0.03`$, for the Ly$`\alpha `$ forest model described in §2.

Remarkably, there should be essentially no pixels in the spectrum at this high redshift where the flux comes close to the continuum level. Only 3% of the spectrum should have a flux higher than 50% of the continuum, and less than 0.5% should have a flux greater than 75% of the continuum \[see the inset of Fig. 2(a)\]. In comparison, Fan et al. (2000) identify the minimum absorption in their observed spectrum to be $`\tau 0.4`$ (although they estimate that continuum uncertainty and noise make this number uncertain by about $`\pm 0.1`$). Absorption of $`\tau 0.4`$ or $`F0.67`$ occurs in approximately 1% of the pixels in our simulated spectra, so we predict about three such pixels to be found in the wavelength range 7300-7956 Å examined by Fan et al. (2000). Therefore, the minimum absorption present in this spectrum is consistent with the predictions of the simple Ly$`\alpha `$ forest model in our simulations. We find that decreasing $`\sigma _8`$ by $`15\%`$ to match the results of McDonald et al. (2000a) increases the inferred value of $`\mathrm{\Gamma }_{12}`$ by $`30\%`$, but does not affect significantly our conclusions about the flux distribution (once the observed $`\overline{F}`$ is matched).

Next, the same lines are plotted in Figure 2(b) at redshift $`z=4.5`$, which is about the mean redshift of the Ly$`\alpha `$ forest in the $`z=5.0`$ quasar spectrum that is shown in Figure 2 of SHCM. The mean transmitted flux in this spectrum is $`\overline{F}=0.25`$, which requires $`\mathrm{\Gamma }_{12}0.21`$ (no error bar was given in SHCM for the value of $`\overline{F}`$). According to SHCM, an upper bound on the minimum optical depth in their spectrum is $`\tau <0.1`$. They identify this region as being produced by a void with density $`\mathrm{\Delta }0.2`$, to derive a lower limit to the background intensity of $`\mathrm{\Gamma }_{12}>2.8`$, for $`\mathrm{\Omega }_bh^2=0.019`$.

Our value of $`\mathrm{\Gamma }_{12}`$, inferred by matching $`\overline{F}`$, is much smaller for the following reasons: Factors of $`(20000/5000)^{0.7}`$, $`(0.4/1.0)^{1/2}`$, and $`0.65/0.5`$ are obtained from our different choices of $`T_0`$, $`\mathrm{\Omega }_m`$, and $`h`$, respectively (note: we assume SHCM used $`h=0.5`$, but their paper is not clear on this point). Our choices place more conservative lower limits on $`\mathrm{\Gamma }_{12}`$, and are more consistent with independent observations. The remaining factor 2.4 must be traced to the identification by SHCM of observed $`\tau =0.1`$ with underdensity $`\mathrm{\Delta }=0.2`$. In our simulations, $`\mathrm{\Delta }<0.2`$ occurs in 1.7% of pixels in real space, and this probability is increased when the expansion of voids in redshift space is considered. Consequently, the minimum absorption in the 200Å region considered by SHCM probably arises from a void with $`\mathrm{\Delta }<0.2`$. Furthermore, the peculiar expansion of voids reduces the optical depth at fixed density. After the inclusion of all these effects through our simulations, and using $`\mathrm{\Gamma }_{12}=0.21`$, we predict $`0.4`$ regions with $`\tau <0.1`$ in the wavelength range $`68007000`$Å identified by SHCM. Lastly, the minimum absorption statistic is biased toward regions with strong noise contribution, so the value $`\tau <0.1`$ is probably not conservative (our simulations give a completely acceptable 1% probability for $`\tau <0.15`$).

Finally, the quasar in Stern et al. (2000) at $`z=5.5`$ has mean transmitted flux $`\overline{F}=0.10\pm 0.02`$. The ionization rate required to match this $`\overline{F}`$ is $`\mathrm{\Gamma }_{12}=0.13\pm 0.03`$, at a mean redshift $`\overline{z}=4.93`$. Their spectrum as shown has only $`15`$ Å resolution so we do not interpret it further.

The evolution of the background intensity $`\mathrm{\Gamma }_{12}`$ derived from the observed mean flux decrement, for $`\mathrm{\Omega }_bh^2=0.02`$, is shown in Figure 3, using the values from the three quasars discussed here and the lower redshift results from McDonald et al. (2000a) (updated by McDonald et al. 2000b using the measured IGM temperatures).

## 4 DISCUSSION

There is a clear reduction of the proper intensity of the ionizing background with redshift, by a factor $`3`$ from $`z=3`$ to $`z=5`$. This reduction was noticed previously by Rauch et al. (1997) up to $`z=4`$. With the new high redshift quasars, we find a continuation of this gradual decrease of $`\mathrm{\Gamma }`$, but without a dramatic drop that would be the sign of approaching the reionization epoch. The consistency of the fraction of the Ly$`\alpha `$ spectrum with low absorption with the observations also suggests that the background intensity fluctuations are not yet extreme at $`z=5`$.

The evolution of the background intensity is related to the evolution of the comoving emissivity, $`ϵ`$, and the proper mean free path $`\lambda `$, by $`ϵ(1+z)^3\mathrm{\Gamma }\lambda ^1`$, if $`\lambda `$ is small compared to the cosmic horizon. Assuming that the rate of incidence of Lyman limit systems evolves as $`(1+z)^{1.5}`$ (Storrie-Lombardi et al. 1994), the mean free path of ionizing photons is $`\lambda (1+z)^4`$, and therefore the comoving emissivity changes between $`z=3`$ and $`z=5`$ by a factor $`(1/3)(6/4)^{43}=0.5`$, indicating only a small reduction, which is clearly different from the reduction by a factor 0.14 predicted by the evolution of the comoving density of luminous quasars (Schmidt, Schneider, & Gunn 1995) and implies, therefore, that galaxies or other sources are dominating the emission of ionizing photons at $`z=5`$.

We thank Nick Gnedin for providing a copy of his Hydro-PM code and for helpful comments on the manuscript, Renyue Cen for the output from his hydrodynamic simulation, and Adam Steed for helpful comments on the manuscript.