text:
# Steps to the Planck Function : A Centenary Reflection

## 1 Introduction

The Planck function was first presented at a meeting of the German Physical Society in Berlin on 19 October 1900 in a contribution entitled ”On an improvement of the Wien spectral equation”, and published in the proceedings of the society in a brief three-page report (Planck, 1900Va). Only the functional form was given, including two undetermined constants, derived on the basis of a tentative relation between the entropy and the energy of each of a collection of radiation resonators. The definitive formula was presented at a meeting of the society on 14 December 1900 in a paper entitled ”On the theory of the law of the energy distribution in the normal spectrum” and published, in a nine-page report (Planck, 1900Vb). Herein Planck introduced his new natural constant $`h`$ relating the ”energy element” $`ϵ`$ associated with a resonator of frequency $`\nu `$ by the formula $`ϵ=h\nu `$, with the numerical value $`h=6.55x10^{27}`$ erg seconds. He only sketches the actual calculation of the entropy, now for the first time related to the probability of resonator complexions via the Boltzmann principle. Thus ended the long search for the blackbody spectral distribution function and simultaneously began the era of quantum physics. A full treatment of the derivation was later given in a paper (Planck, 1901a) entitled ”On the law of the energy distribution in the normal spectrum” submitted to Annalen der Physik on 9 January 1901.

## 2 Thermal Radiation

That the arrival of both the visible light and the ’dark’ radiant heat from the Sun are simultaneously interrupted and restored at the time of a solar eclipse, indicated that light and heat are propagated across space with the same velocity. The astronomer W. Herschel (1800), using a thermometer exposed to the prismatic spectrum of the Sun, discovered that the heating properties of light extended beyond the visible spectrum, being enhanced most especially in the infrared. He also showed that the unseen heat radiation obeyed the same laws of reflection and refraction as visible light. The similarity in other properties, including bi-refringence and polarisation, as found by others, for light and heat, illustrated a difference of degree rather than kind. Furthermore, that sources of light are also sources of heat, and $`\mathrm{v}\mathrm{i}\mathrm{c}\mathrm{e}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{a}`$, served to distinguish the two manifestations only by the means of detection employed. Indeed, from his measurements, Herschel presented two distributions of the (continuous) intensity of the solar rays, based respectively on their ability to illuminate and to heat. Because of the two different techniques of measurement, the relative scales are not related. If the two intensity distributions were brought to a common scale and then summed, it would represent the first wideband total intensity distribution for the solar spectrum.

Simple experiments, such as those by Leslie (1804) and others, showed that the rate of emission of thermal radiation from a heated body, expressed as the emissive power or radiant emittance E, being the power emitted per unit area, depends on the temperature and on the nature of the surface. Thus a blackened surface has a higher emittance than a polished one at the same temperature. Also when two such surfaces are subjected to incident radiant energy the fraction of energy absorbed, or absorptivity $`A`$, depends on the temperature and on the nature of the surface, a blackened surface being a more efficient absorber than a polished one.

Kirchhoff (1859), by means of a thought experiment, demonstrated that the ratio $`E/A`$ is the same for all bodies at the same temperature, emitting and absorbing at the same wavelength. He considered two bodies, in the form of infinte surfaces, S1 and S2, backed with perfectly reflecting mirrors and arranged parallel to each other with the mirrors on the outside. S1 is such that it can only emit and absorb radiation of wavelength $`\mathrm{\Lambda }`$, while S2 can emit and absorb all possible wavelengths $`\lambda `$. At wavelength $`\mathrm{\Lambda }`$ the emissivities are $`E_1`$ and $`E_2`$, and the absorptivities are $`A_1`$ and $`A_2`$. For all $`\lambda \mathrm{\Lambda }`$ all the radiation emitted by S2 is reabsorbed by it. Therefore for temperature equilibrium S2 must absorb as much as it emits of the radiation of wavelength $`\lambda =\mathrm{\Lambda }`$, whether it comes from S1 or S2. Of the radiation $`E_1`$ emitted by S1, S2 absorbs, after repeated absorptions and emissions :

$$A_2E_1\underset{n=0}{\overset{\mathrm{}}{}}k^n=A_2E_1/(1k)$$

while of the radiation $`E_2`$ emitted by S2 itself, S2 absorbs:

$$(1A_1)A_2E_2\underset{n=0}{\overset{\mathrm{}}{}}k^n=(1A_1)A_2E_2/(1k)$$

where $`k=(1A_1)(1A_2)`$. Hence :

$$E_2=A_2[E_1+(1A_1)E_2]/(1k)$$

giving :

$$A_1E_2=A_2E_1\mathrm{or}E_1/A_1=E_2/A_2$$

The same result is obtained on applying the condition of temperature equilibrium to S1, yielding :

$$E_1=A_1[E_2+(1A_2)E_1]/(1k)$$

Since the composition of the surfaces, as well as the wavelength $`\mathrm{\Lambda }`$ and the temperature are arbitrary, it follows that the ratio $`E/A`$, at any wavelength or temperature, is the same for all surfaces, and is a function only of wavelength and temperature, which statement is known as Kirchhoff’s Law.

Kirchhoff (1860) generalised the above argument, which involved only discrete wavelengths, to the continuous spectrum, by re-defining $`E`$ as the emittance per unit wavelength interval. He also introduced the concept of a perfectly black body, being one which absorbs all incident radiation, for which $`A=1`$ at all wavelengths and temperatures. By considering the temperature equilibrium of a black body enclosed in a cavity with walls made of a similarly black material, he showed that the emittance is equal to the radiant power, $`F(\lambda ,T)`$ per unit wavelength interval, crossing unit area, within the ‘black body’ enclosure. Hence, by considering the temperature equilibrium of a non-black body within the cavity, Kirchhoff’s Law takes the form $`E/A=F(\lambda ,T)`$. Thus the absorptivity of a body is equal to its relative emittance, or the emittance of the body relative to that of a black body. Kirchhoff stressed the importance of the study of black body radiation in order to determine the function $`F(\lambda ,T)`$. It is also worth noting here that the enunciation by Kirchhoff of the principle of detailed balance between emission and absorption precedes that of Einstein by more than half a century.

## 3 Integrated Properties of Black Body Radiation

Given the dependence on wavelength $`\lambda `$ and temperature $`T`$ of the black body radiation intensity $`I(\lambda ,T)`$, its energy density $`u`$ and flux $`F`$, both monochromatic and integrated, can be obtained :

$$I(T)=I(\lambda ,T)d\lambda $$

$$u(\lambda ,T)=(4\pi /c)I(\lambda ,T),u(T)=u(\lambda ,T)d\lambda =(4\pi /c)I(T)$$

$$F(\lambda ,T)=\pi I(\lambda ,T),F(T)=F(\lambda ,T)d\lambda =(c/4)u(T)$$

As indicated by Kirchhoff, the properties of black body radiation can be studied by observations on the radiation emerging from a small orifice in the wall of a black body cavity. Experiments by Stefan (1879) established that the total (integrated) emergent flux could be represented in the form :

$$F(T)=\sigma T^4$$

where $`\sigma =5.670x10^8Jm^2s^1K^4`$ is the Stefan(-Boltzmann) constant, from which it follows that the energy density is given by :

$$u(T)=(4\sigma /c)T^4=aT^4$$

with the radiation constant $`a=7.565x10^{16}Jm^3K^4`$. The latter form was established theoretically by Boltzmann (1884) on the basis of thermodynamics. Using for the pressure of radiation the relation :

$$P(T)=(1/3)u(T)$$

as indicated by the electromagnetic theory of Maxwell (1873), and assuming that radiation obeyed the known laws of thermodynamics, then putting $`U=uV`$ :

$$\mathrm{d}Q=T\mathrm{d}S=\mathrm{d}(uV)+P\mathrm{d}V=(u+P)\mathrm{d}V+V\mathrm{d}u$$

from which it follows that :

$$(/u)_V[(u+P)/T]=(/V)_u[V/T]$$

which, on using $`P=u/3`$, gives :

$$(T/u)(u/T)_V=(\mathrm{ln}u/\mathrm{ln}T)_V=4$$

so that :

$$u=aT^4,P=(1/3)aT^4$$

in agreement with the result of Stefan. Alternatively, Stefan’s result may be taken as a proof of the relation $`P=u/3`$. From the above expressions for $`u`$ and $`P`$ it follows that :

$$\mathrm{d}S=(4/3)aT^3\mathrm{d}V+4aT^2V\mathrm{d}T$$

giving for the entropy :

$$S=(4/3)aT^3V$$

## 4 Spectrum of Black Body Radiation

The challenge, first presented by Kirchhoff (1859), to determine the spectral distribution of black body radiation, was not to be completely answered until the turn of the century. Experiments by Draper (1847) and by Weber (1887) showed that, with increasing temperature, the maximum intensity in the light emitted by incandescent solid bodies, progressed across the spectrum, in the direction from the red to the violet. Bartoli (1876) had pointed out that, in order not to contradict the second law of thermodynamics, radiation must exert a pressure, and Boltzmann (1884) showed that work must be done against this pressure in compressing a volume occupied by radiation alone. Wien (1893), developing this argument, made the first important steps in the determination of the black body radiation spectrum. On the basis of the Doppler principle, applied to radiation within a rectangular parallelepiped and reflected (normally) from a wall moving with velocity $`v`$, where the relative wavelength change is given by $`\delta \lambda /\lambda =2v/c`$ (twice the ordinary Doppler shift), then radiation of wavelength $`\lambda _0`$ and frequency $`\nu _0`$ in a volume $`V_0`$, subjected to an adiabatic compression to volume $`V`$, undergoes a change to wavelength $`\lambda `$ and frequency $`\nu `$ where :

$$\lambda ^3V^1=\lambda _0^3V_0^1\mathrm{and}\nu ^3V=\nu _0^3V_0$$

The change in energy of the radiation, supplied by the work done in the compression against the radiation pressure, results in a change in temperature from $`T_0`$ to $`T`$. Since for an adiabatic change $`VT^3=V_0T_0^3`$ it follows that :

$$\lambda T=\lambda _0T_0\mathrm{and}\nu /T=\nu _0/T_0$$

giving the earliest expression of the Wien Displacement Law, which states that with changing temperature the black body spectrum is shifted so that for each wavelength the product of the temperature and the wavelength remains constant. Similarly each interval of wavelength and frequency is changed according to :

$$T\mathrm{d}\lambda =T_0\mathrm{d}\lambda _0\mathrm{and}\mathrm{d}\nu /T=\mathrm{d}\nu _0/T_0$$

Since the (spatial) energy density obeys the relation :

$$u(\lambda ,T)\mathrm{d}\lambda /u(\lambda _0,T_0)\mathrm{d}\lambda _0=u(\nu ,T)\mathrm{d}\nu /u(\nu _0,T_0)\mathrm{d}\nu _0=T^4/T_0^4$$

then :

$$u(\lambda ,T)/u(\lambda _0,T_0)=(T/T_0)^5=(\lambda _0/\lambda )^5$$

$$u(\nu ,T)/u(\nu _0,T_0)=(T/T_0)^3=(\nu /\nu _0)^3$$

in agreement with the observations of Weber. It also follows that the wavelength $`\lambda _m(T)`$ and frequency $`\nu _m(T)`$ at which the wavelength and frequency spectral distributions reach their maxima satisfy the relations :

$$T\lambda _m(T)=T_0\lambda _m(T_0)\mathrm{and}\nu _m(T)/T=\nu _m(T_0)/T_0$$

representing another expression of the Wien Displacement Law. Therefore given the spectral distribution for one temperature $`T_0`$, that for any other temperature $`T`$ can be determined.

Wien (1894) re-derived all the above relations when the radiation is enclosed in a spherical volume of variable radius. He also showed that the spectral distributions for different temperatures cannot cross. The condition for this to be true is that :

$$[\mathrm{ln}u(\lambda ,T)/\mathrm{ln}\lambda ]_T>5$$

which is trivially true for the rising branch of the curve, and also obviously true for the falling branch. In terms of the frequency distribution the condition becomes :

$$[\mathrm{ln}u(\nu ,T)/\mathrm{ln}\nu ]_T<+3$$

### 4.1 The Wien Formula

Wien (1896) pointed out that, since the black body radiation spectrum corresponded to thermal equilibrium, it should be completely determined from the principle of maximum entropy, although he did not do this himself. Instead, adopting a radiation model based on analogy with the Maxwell velocity distribution, he obtained a spectrum of the form :

$$u(\lambda ,T)=f_1(\lambda )\mathrm{exp}[f_2(\lambda )/T]$$

where $`f_1(\lambda )=c_1/\lambda ^\alpha `$, ($`c_1`$ and $`\alpha `$ being constants), and $`f_2(\lambda )`$ are functions to be determined. In order to satisfy the Stefan-Boltzmann $`T^4`$ dependence of the energy density, it is necessary for $`\alpha =5`$ and for $`f_2(\lambda )/T`$ to be a function of the product $`\lambda T`$, thus giving a more general mathematical expression to the Wien Displacement Law. Putting, specifically, $`f_2(\lambda )/T=c_2/(\lambda T)`$, where $`c_2`$ is a constant, gives the form :

$$u(\lambda ,T)=c_1\lambda ^5\mathrm{exp}[c_2/(\lambda T)]$$

known as Wien’s formula, which attains its maximum value at a wavelength $`\lambda _m`$ with $`\lambda _mT=c_2/5`$. The formula was in good accord with the observations of Paschen (1896), with $`c_2=5\lambda _mT=14455\mu `$K, or $`\lambda _mT=2.891x10^3`$ m K.

In a series of studies, from 1897 onwards, on the entropy of radiation, treating the radiation field as an ensemble of monochromatic resonators, Planck sought to realize the hope of Wien that the blackbody radiation spectrum could be derived on the basis of the maximization of the entropy. Indeed Planck (1899) believed (erroneously) that he had shown that the Wien spectral distribution was a necessary consequence of the principle of entropy maximization, and that its range of validity coincided with that of the (two) laws of thermodynamics. Such support for the Wien formula led to it being referred to as the Wien-Planck distribution. However, experimental measurements by Lummer and Pringsheim (1899), Rubens (1899) and others pointed to systematic departures of the Wien-Planck formula from the experimental results particularly at long wavelengths. Thiesen (1900), following the lead of Wien, expressed the spectral distribution in the, only slightly different, general form :

$$u(\lambda ,T)=T^5\psi (\lambda T)$$

where $`\psi (\lambda T)`$ is a still unknown function. Thiesen also suggested that the function $`\psi `$, when determined, would accommodate the deviations between theory and experiment by differing from the Wien-Planck form, leading to a modification of the Wien displacement law.

Planck advanced his support of the Wien formula in two further papers (Planck, 1900Aa and 1900Ab). In the first paper he adopted, by way of definition, $`s=(w/a\nu )\mathrm{ln}(w/eb\nu )`$, where $`s`$ and $`w`$ are the entropy and (average) energy per radiation resonator of frequency $`\nu `$, and $`a`$ and $`b`$ are universal positive constants, while $`e`$ is the base of natural (Napierian) logarithms. Using the thermodynamic relation $`ds/dw=1/T`$ where $`T`$ is the temperature, it follows that $`w=b\nu \mathrm{exp}(a\nu /T)`$, from which, on multiplying by the density $`8\pi \nu ^2/c^3`$ of resonators in $`\nu `$-space, the Wien formula follows. From the measurements of Paschen (1899) the constants $`a`$ and $`b`$ were determined to be $`a=4.818x10^{11}`$ sec K and $`b=6.885x1027`$ erg sec. The value of $`b`$ here therefore represents the first determination of what later became the Planck constant $`h`$, while the value of $`a`$ is that of $`h/k`$. In the second paper Planck sought to justify his earlier assumed definition of the entropy $`s`$, by arguing that maximization required $`d^2s/dw^2=F(w)`$, where $`F(w)`$ is a positive function of $`w`$. Using a further relation (of limited validity) $`F(nw)=F(w)/n`$, where $`n`$ is the number of resonators, Planck adopted the solution $`F(w)1/w`$. Thence putting $`d^2s/dw^2=\alpha /w`$, it follows that $`s=\alpha w\mathrm{ln}(\beta w/e)`$, where now $`\alpha `$ and $`\beta `$ are positive constants which depend only on the frequency $`\nu `$. Then $`w=[1/e\beta (\nu )]\mathrm{exp}[1/T\alpha (\nu )]`$, leading again to a spectral distribution of the Wien form. Furthermore the Wien displacement law required both $`\alpha `$ and $`\beta `$ to be proportional to $`\nu ^1`$, thus permitting the identifications $`\alpha =1/a\nu `$ and $`\beta =1/eb\nu `$, where $`a`$ and $`b`$ are the quantities introduced earlier. However further work by Lummer and Pringsheim (1900) showed that both $`a`$ and $`b`$ were not constant but increased with wavelength, leading them to conclude the invalidity of the Wien-Planck spectral distribution.

### 4.2 The Rayleigh Formula

To Rayleigh (1900) the Wien formula was difficult to accept on the grounds it implied that, with increasing temperature, the radiation density at a given wavelength approached a limit. Using an argument based on the still disputed Maxwell-Boltzmann principle of the equipartition of energy, in which each degree of freedom has an energy proportional to the temperature $`T`$, and in analogy with the theory of sound, taking for the density of vibrational modes $`n(\lambda )\lambda ^4`$, he obtained :

$$u(\lambda ,T)=KT\lambda ^4$$

where $`K`$ is a constant, which has the form prescribed by Thiesen, suggesting that it may be correct for large values of $`\lambda T`$. However the formula tends to infinity at short wavelengths, leading Rayleigh to attach the Wien exponential factor $`\mathrm{exp}[c_2/(\lambda T)]`$.

### 4.3 The Planck Formula

In the light of the increasing evidence for the inadequacy of the Wien blackbody spectral distribution, Planck (1900Va) now addressed the problem of improving upon it. He recalled his earlier conclusion that the maximization of entropy was not of itself sufficient to derive the correct distribution, and recognised that his support of the Wien function had depended upon his adoption of the further special condition $`d^2s/dw^2=\alpha /w`$, applicable only in the case of indefinitely small values of $`w`$ and $`s`$. Accordingly he considered instead the condition $`d^2s/dw^2=\alpha /w(\beta +w)`$, which he described as by far the next simplest relation giving $`s`$ as a logarithmic function of $`w`$, while at the same time reducing to the former relation for vanishingly small values of $`w`$. It is readily shown that, ignoring constants of integration, $`ds/dw=(\alpha /\beta )\mathrm{ln}\{w/(\beta +w)\}`$, and $`s=(\alpha /\beta )\{w\mathrm{ln}w(\beta +w)\mathrm{ln}(\beta +w)\}`$. Then, using $`ds/dw=1/T`$, it follows that $`w=\beta /\{\mathrm{exp}(\beta /\alpha T)1\}`$. Noting that the Wien Displacement Law requires (see later) that $`s`$ be a function of $`w/\nu `$, where $`\nu `$ is the frequency, (demanding that $`\beta `$ be proportional to $`\nu `$), then enabled Planck to express the spectral distribution in the form $`u(\lambda ,T)=C_1\lambda ^5/\{\mathrm{exp}(C_2/\lambda T)1\}`$, where $`C_1`$ and $`C_2`$ are constants, thus introducing what is now known as the Planck Function. Rubens and Kurlbaum (1900) pointed out that the new Planck formula included those of Wien and Rayleigh in the limiting cases of short and long wavelengths respectively, while their experimental measurements favoured it over all others.

However, the derivation was largely phenomenological since the entropy had not been calculated from first principles. Therefore Planck (1900Vb) next sought to follow the ”Boltzmann prescription” whereby entropy involves disorder through the relation $`S=k\mathrm{ln}W`$, where $`W`$ is the number of ways or states in which the physical system may be realized. This use of this famous relation, never given by Boltzmann himself yet inscribed on his tombstone, is the first occurrence known to the writer. Since the total radiation energy is distributed among all the radiation resonators involving a number of frequencies, it is first necessary to consider the distribution of an energy $`E`$ among the $`N`$ resonators with the same frequency $`\nu `$. Now if $`E`$ were infinitely divisible the distribution would be possible in infinitely many ways. Hence Planck considered \- and this is the essential point of the whole calculation - $`E`$ to be made up of a fixed number $`P`$ of finite equal parts, or energy elements, each of magnitude $`ϵ=h\nu `$, where $`h`$ is a new constant of nature, so that $`E=Pϵ=Nw`$, where $`w`$ is the average energy per resonator. This therefore represents the first appearance of the Planck constant, both as to its symbol and definition. Planck also proceeded to use combination theory to calculate the total number of possible complexions :

$$W=(N+P1)!/(N1)!/P!(N+P)^{N+P}/N^N/P^P$$

on making use of the Stirling formula whereby for large $`n,n!n^n`$. He did not however give an expression for the entropy $`S`$, but gave his final result for the spectral distribution in its definitive form :

$$u(\nu ,T)=(8\pi h\nu ^3/c^3)/[\mathrm{exp}(h\nu /kT)1]$$

where now the Boltzmann constant $`k`$ appears explicitly as well as the new (Planck) constant $`h`$. From the formula, using the measurements of Kurlbaum (1898) and of Lummer and Pringsheim (1900), the values of $`h`$ and $`k`$ were determined to be $`h=6.55x10^{27}`$ erg sec and $`k=1.346x10^{16}`$ erg grad.

Very shortly afterwords, Planck (1901a) completed his derivation with details of the calculation of the entropy giving :

$$\mathrm{ln}W=N[(1+P/N)\mathrm{ln}(1+P/N)(P/N)\mathrm{ln}(P/N)]$$

so that, since $`P/N=w/ϵ`$, the entropy per resonator becomes :

$$s=S/N=k[(1+w/ϵ)\mathrm{ln}(1+w/ϵ)(w/ϵ)\mathrm{ln}(w/ϵ)]=f(w/ϵ)$$

which for small values of $`w/ϵ`$ (small $`w`$ and large $`ϵ`$ or $`\nu `$) leads to the earlier (Wien) result. Now it follows that :

$$\mathrm{d}s/\mathrm{d}w=1/T=(k/ϵ)\mathrm{ln}(1+ϵ/w)$$

giving :

$$w=ϵ/\{\mathrm{exp}(ϵ/kT)1\}$$

Then since :

$$u(\nu ,T)=8\pi (\nu ^2/c^3)w$$

and, according to the Wien Displacement Law :

$$u(\nu ,T)=\nu ^3f_\nu (\nu /T)$$

where $`f_\nu `$ is a still undetermined function :

$$w=\nu [(c^3/8\pi )f_\nu (\nu /T)]$$

Thus :

$$1/T=\mathrm{d}s/\mathrm{d}w=\nu ^1(8\pi /c^3)f_\nu ^1(w/\nu )$$

so that, on integrating :

$$s=F(w/\nu )$$

and therefore, in order to agree with the result that $`s=f(w/ϵ)`$, it follows that $`ϵ=h\nu `$ giving :

$$w=h\nu /[\mathrm{exp}(h\nu /kT)1]=(hc/\lambda )/[\mathrm{exp}(hc/kT\lambda )1]$$

where $`h`$ is again the Planck constant. Thus one finally obtains the Planck formula :

$$u(\nu ,T)=(8\pi h\nu ^3/c^3)/[\mathrm{exp}(h\nu /kT)1]$$

$$u(\lambda ,T)=(8\pi hc/\lambda ^5)/[\mathrm{exp}(hc/kT\lambda )1]$$

The integrated radiation energy density is now $`u(T)=aT^4`$ with (Planck, 1901a) :

$$a=(8\pi k^4/h^3c^3)_0^{\mathrm{}}z^3/(\mathrm{exp}z1)dz=8\pi ^5k^4/15c^3h^3$$

in which the integral evaluates to $`3!\zeta (4)=\pi ^4/15`$ giving $`a=7.565x10^{16}m^3K^4`$. From the expression for w, giving $`s=k[z/(\mathrm{exp}z1)+\mathrm{ln}\{\mathrm{exp}z/(\mathrm{exp}z1)\}]`$, where $`z=h\nu /kT`$, the integrated entropy density is $`S(T)=(4/3)aT^3V`$, as has already been derived from thermodynamics.

The maximum of the wavelength distribution occurs for $`x_m=hc/kT\lambda _m`$ where :

$$x_m=5[1\mathrm{exp}(x_m)]=4.965114\mathrm{}$$

while the maximum of the frequency distribution occurs for $`y_m=h\nu _m/kT`$ where :

$$y_m=3[1\mathrm{exp}(y_m)]=2.821439\mathrm{}$$

giving $`\nu _m\lambda _m/c=y_m/x_m=0.568252\mathrm{}`$, indicating that $`\lambda _m<c/\nu _m`$ and $`\nu _m<c/\lambda _m`$. Hence also :

$$\nu _m/T=(k/h)x_m=5.8785x10^{10}HzK^1$$

$$\lambda _mT=(hc/k)/y_m=2.89779x10^3mK$$

In the low frequency/long wavelength limit $`h\nu /kT=hc/kT\lambda 1`$ the Rayleigh formula is obtained, on expanding the exponential in the denominator :

$$u(\nu ,T)=(8\pi h\nu ^2/c^3)kT$$

$$u(\lambda ,T)=8\pi kT/\lambda ^4$$

while in the high frequency/short wavelength limit $`h\nu /kT=hc/kT\lambda 1`$ the Wien formula is obtained, on ignoring the 1 in the denominator :

$$u(\nu ,T)=(8\pi h\nu ^3/c^3)\mathrm{exp}(h\nu /kT)$$

$$u(\lambda ,T)=(8\pi hc/\lambda ^5)\mathrm{exp}(hc/kT\lambda )$$

Further experimental measurements by Lummer and Pringsheim (1901) showed that the values of $`\lambda _mT`$ were in excellent agreement with the predictions of Planck. Finally, Planck (1901b) returned to the question of the entropy of radiation, to demonstrate that the expression obtained from his theory satisfied the maximization principle. Thus in his stumbling search for the blackbody radiation spectral distribution, Planck, in forging the hitherto missing link between radiation and thermodynamics, not only closed the era of classical physics but also opened up the era of quantum physics.

## 5 Biographical Note

Max (Karl Ernst Ludwig) Planck was born on 23 April 1858 in Kiel, and died on 4 October 1947 in Göttingen. He attended the universities of Munich and Berlin (in the latter he studied under Kirchhoff and Helmholtz), and received his doctorate from the University of Munich in 1879. Subsequently he taught in the universities of Munich (1880-1885), Kiel (1885-1889) and Berlin (1889-1928). In 1894 he became a member of the Prussian Academy of Science, and the permanent secretary of its Section of Mathematics and Natural Science in 1912. He was awarded the Nobel Prize for Physics in 1918. In 1930 he became president of the Kaiser Wilhelm Society in Berlin, later to be known as the Max Planck Society.

## 6 References

Bartoli, X., 1876, Fortschritte (2) 32, 888 and 1541.
Boltzmann, L., 1884, Annalen der Physik 22, 31.
Draper, J.W., 1847, Phil. Mag. 30, 345.
Herschel, W., 1800, Phil. Trans. Roy. Soc. London, pages 255, 284, 293, 437 (also 1800, Phil. Mag. 7, 311 and idem 8, pages 9, 16, 126, 253).
Kirchhoff, R., 1859, Monatsbericht d. Berl. Akad., December 1859, 783.
Kirchhoff, R., 1860, Annalen der Physik 109, 275.
Kurlbaum, F., 1898, Wied. Ann. 65, 759.
Leslie, J., 1804, Experimental Inquiry into the Nature and Propagation of Heat.
Lummer, O. & Pringsheim, E., 1899, Verh. d. deut. phys. Ges. 1, 23 and 215.
Lummer, O. & Pringsheim, E., 1900, Verh. d. deut. phys. Ges. 2, 163.
Lummer, O. & Pringsheim, E., 1901, Verh. d. deut. phys. Ges. 3, 36.
Maxwell, J. C., 1873, A Treatise on Electricity and Magnetism.
Paschen, F., 1899, Sitz. d. k. Akad. d. Wiss. zu Berlin, 405.
Planck, M., 1900Aa, Annalen der Physik 1, 69.
Planck, M., 1900Ab, Annalen der Physik 1, 719.
Planck, M., 1900Va, Verh. d. deut. phys. Ges. 2, 202.
Planck, M., 1900Vb, Verh. d. deut. phys. Ges. 2, 237.
Planck, M., 1901a, Annalen der Physik 4, 553.
Planck, M., 1901b, Annalen der Physik 6, 818.
Rayleigh, Lord (Strutt, J. W.), 1900, Phil. Mag. 49, 539.
Rubens, H., 1899, Wied. Ann. 69, 582.
Rubens, H. & Kurlbaum, F., 1900, Sitz. d. k. Akad. d. Wiss. zu Berlin, 929.
Stefan, J.,1879, Sitz. d. k. Akad. d. Wiss. zu Wien, 79 (II), 391.
Thiesen, M., 1900, Verh. d. deut. phys. Ges. 2, 65.
Weber, H. F., 1887, Sitz. d. k. preuss. Akad. d. Wiss. zu Berlin, 491.
Wien, W., 1893, Sitz. d. k. preuss. Akad. zu Berlin, 55.
Wien, W., 1894, Annalen der Physik 52, 132.
Wien, W., 1896, Annalen der Physik 58, 662 (also 1897, Phil. Mag. 43, 214).


format: text

text:
# Stellar Pollution in the Solar Neighborhood

## 1 INTRODUCTION

The discovery that at least $`68\%`$ of solar type stars harbor Jupiter-mass or larger bodies, often in small, eccentric orbits (Mayor & Queloz, 1995; Marcy & Butler, 1996; Butler & Marcy, 1996; Marcy et al., 2000) shows that planets are not exceptionally rare. The transiting planet orbiting HD 209458 has a mass of $`0.69`$ Jupiter masses, and a radius about 1.4 times that of Jupiter. Clearly it is a gas giant like Jupiter or Saturn; the minimum masses of the other known objects suggest that they are also gas giants. The Doppler technique used for most of the discoveries cannot find terrestrial mass objects in AU scale orbits. What fraction of solar type stars have terrestrial-type planets?

How could an observer with our current technology located on a star in the solar neighborhood decide if there were any terrestrial-type bodies orbiting the sun? More generally, how could such an observer estimate the fraction of solar type stars having companions of terrestrial (as opposed to gas giant) compositions? In this paper we explore one possible way of addressing the latter question, stellar-mass dependent photospheric metallicity trends.

In section 2 we argue that a few Earth masses of rocky/icy material have accreted onto the Sun over its lifetime. In section 3 we discuss one possible observational tracer of similar accretion occurring on other solar type stars. In section 4 we examine a sample of stars with spectroscopically measured photospheric iron abundances

$$[\mathrm{Fe}/\mathrm{H}]\mathrm{log}\left[\frac{f_{Fe}}{f_{Fe,}}\right],$$
(1)
where $`f_{Fe}`$ is the mass abundance of iron in the photosphere of the star, and $`f_{Fe,}1.3\times 10^3`$ is the the mass abundance of iron in the photosphere of the sun. We compare the observations to Monte Carlo models of stellar pollution in section 5. We discuss the implications of our results in section 6 and present our conclusion in the final section.

## 2 POLLUTION BY SMALL BODIES IN OUR SOLAR SYSTEM

The surface density $`\mathrm{\Sigma }_{Fe}`$ of iron in our solar system follows a rough power law with distance from the sun Weidenschilling (1977). We have compiled recent estimates for the iron content of the planets in table 1. From a least squares fit including Venus, Earth, and the four gas giants we find $`\mathrm{\Sigma }_{Fe}4(r/1AU)^{1.7}\mathrm{g}\mathrm{cm}^2`$. Weidenschilling noted a clear deficit of material in the region interior to Venus, and between Earth and Jupiter, relative to this power law. Gas drag acting on planetesimals near the sun could reduce the surface density of such bodies near the present orbit of Mercury by dropping them onto the sun. Because the solar convection zone was very deep for $`20`$ million years (see figure 1), while the gas disk that would have produced the drag is believed to have survived only a few million, this would not have altered the apparent metallicity of the sun.

Materiel in the asteroid belt clearly survived such gas drag. However, Kirkwood noticed over one hundred years ago that the distribution of asteroids showed distinct gaps at the location of orbital resonances with Jupiter Kirkwood (1867). This suggests that material has been removed from the gaps under the influence of Jupiter. Over the last twenty years, the dynamics of this removal have been worked out in considerable detail Wisdom (1983); Holman & Murray (1996); Gladman et al. (1997). The gaps are the result of resonant, chaotic perturbations of the asteroid orbits by Jupiter. A second feature, seen in plots of orbital eccentricity versus semimajor axis, is the result of a secular resonance, the $`\nu _6`$ resonance. In either type of resonance bodies with semimajor axis larger than about half that of Jupiter’s tend to be removed from the solar system; their eccentricities grow with time, increasing their apoapses until they reach the orbit of Jupiter. Close encounters with that planet then quickly eject the asteroids.

Resonant asteroids in smaller orbits tend not to be ejected; their apoapse does not reach the orbit of Jupiter. Instead, their periapses decrease until they hit the sun. Several related mechanisms can combine to produce this decrease in periapse, but the two main actors are mean motion resonances and the $`\nu _6`$ resonance already mentioned. In a mean motion resonance such as that responsible for the $`3/1`$ Kirkwood gap, the eccentricity of the resonant asteroid undergoes a random walk. Since the eccentricity cannot decrease below zero, there is a tendency for bodies with small $`e`$ to drift to larger $`e`$. When the eccentricity is large enough, the asteroid strikes the sun Gladman et al. (1997). A small fraction will suffer close encounters with the Earth or other terrestrial bodies; this may lead to Jupiter crossing orbits, which usually leads to ejection of the asteroid from the solar system. Only a very small fraction will collide with either a terrestrial planet or Jupiter. In the $`\nu _6`$ secular resonance the precession rate equals the sixth secular frequency of the solar system. This frequency appears with substantial amplitude in the orbital motion of both Jupiter and Saturn; it is roughly the average precession rate of the latter planet. The $`\nu _6`$ resonance rapidly removes angular momentum from the orbit of the asteroid, dropping the periapse below the surface of the sun.

These resonances are currently active, and are believed to be responsible for the bulk of the meteorites striking the Earth. While there is no direct observational confirmation, it is virtually certain that a much larger flux of asteroidal material from these resonances is currently hitting the sun.

It is also likely that the flux of asteroidal material was much larger in the past. The locations of both mean motion and secular resonances depends on the positions of Jupiter and Saturn, and on any massive gas disk. It has been suggested that as the protoplanetary gas disk dissipated the $`\nu _6`$ resonance swept through the belt, depleting much of the mass. More recently it has been suggested that one or both of Jupiter and Saturn migrated to their current positions, causing a similar sweeping of mean motion resonances. This would preferentially deplete the outer belt Holman & Murray (1996); Liou & Malhotra (1997), but would also drop material from the inner belt onto the sun. Finally, we note that the asteroid belt is hot, in the sense that the distribution of both $`e`$ and $`i`$ is broad. A number of authors have suggested that a planet size body swept through the region early in the history of the solar system, heating the belt as well as dropping material on the star Wetherill (1992); Petit et al. (1991).

The surfaces of the terrestrial planets and the Moon show evidence for an extended period of bombardment, known as the late heavy bombardment, lasting up to a billion years after their surfaces cooled. This is consistent with some of the longer lived depletion scenarios described above.

From table 1, the amount of material inferred to have populated the original asteroid belt is of order five to ten Earth masses. Roughly 2-5 Earth masses would have been between the present position of Mars and $`2.5`$ AU, half the present semimajor axis of Jupiter. A substantial fraction of this material, of order half, would have ended up in the sun. We will take $`2M_{}`$ as a representative number. Meteoritic material is roughly $`20\%`$ iron by weight Grevesse & Anders (1989), so we take $`0.5M_{}`$ as the amount of iron that was dropped on the sun after the convection zone thinned.

## 3 POSSIBLE OBSERVATIONAL EVIDENCE FOR POLLUTION

There is another possible piece of evidence besides the late heavy bombardment for such a late depletion of the asteroid belt. As we noted above, a substantial fraction of the “missing” material originally in the asteroid belt strikes the sun. This material will be mixed throughout the convection zone of the star. Since a substantial bombardment of the terrestrial planets lasted for much longer than 100 million years, this material would have landed in a convection zone of mass $`3\times 10^2M_{}`$. Hence there is a possibility that the signature of solar bombardment was left in the form of an enhanced metallicity in the solar envelope.

The mass of the convection zone was $`3\%`$ of the solar mass between $`10^8`$ and $`10^9`$ years, while the present mass fraction of iron observed in the photosphere is $`1.3\times 10^3`$. This yields a total iron content of $`12.75M_{}`$ between $`10^8`$ and $`10^9`$ years. Adding half an Earth mass of iron would increase the observed $`[\mathrm{Fe}/\mathrm{H}]`$ by only $`0.017`$ dex.

Precision measurements of the bulk metallicity of a single solar type star are not presently available. Measurements of the sound speed in the interior of the sun, using the five minute oscillations, could in principle provide such a measurement. The only work along these lines that we are aware of is the paper by Henney & Ulrich (1998). They show that the addition of $`8M_{}`$ of cometary material results in a frequency shift that, because of uncertainties in solar models, is too small to be detected.

However, it is possible that observations of large numbers of stars could reveal the presence of pollution. For example, stars having a fixed metallicity but different masses will mix material dropped onto their surfaces to different depths. Knowledge of the depth of this surface mixing layer, which we summarize in the next subsection, allows us to look for variations of metallicity with stellar mass produced by accreted material.

### 3.1 The Mass Of The Surface Mixing Region

Material dropped onto the surface of a star will not in general be confined to the photosphere. The surface layers of stars less than about $`1.4M_{}`$ are convectively unstable, and any material dropped on the star is expected be mixed throughout the convection zone. Observations of lithium (Li) abundance provide support for this theoretical picture. But these same observations suggest that there is another form of mixing in stars more massive than $`1.2M_{}`$, possibly related to meridional circulation.

The evidence for the two types of mixing comes in the form of variations in the surface lithium abundance with stellar mass. Lithium is destroyed when it is exposed to temperatures above about $`3.1\times 10^6K`$; the temperature at the photosphere is much lower than this, so lithium depletion indicates that photospheric material is mixed down to regions where the temperature is higher. Lithium could also be removed from the photosphere by simple settling (so that the lithium is stored below the photosphere, as opposed to being destroyed) or by stellar mass loss. However, the observations, which we review below, point to destruction by mixing.

Observations show that for stars with mass below $`1.2M_{}`$ the photospheric abundance of lithium decreases rapidly with decreasing stellar mass, e.g., Boesgaard (1991). Some depletion occurs before the main sequence, but there is strong evidence that depletion also occurs on the main sequence Jones et al. (1999). The main sequence depletion is believed to be related to the increasing depth of the bottom of the surface convection zone with decreasing stellar mass in these stars, supplemented by some convective overshoot and/or settling in stars near the upper end of this mass range.

Stars above about $`1.2M_{}`$ have very thin, or even lack, convection zones. In standard stellar models, with no settling, mixing or mass loss, the surface Li and Be abundances are constant over the main sequence lifetime of such high mass stars. Both Li and Be are destroyed below modest depths in the interiors of such stars. For stars in the mass range $`1.2`$-$`1.6M_{}`$ and having solar metallicity, our stellar models indicate that the nominally undepleted region has a mass of $`3\times 10^2M_{}`$ for Li and $`6\times 10^2M_{}`$ for Be. (We assume that the lithium is destroyed at temperatures above $`3.1\times 10^6`$K). At depths below those corresponding to these masses, Li and Be should be destroyed.

Real stars with masses larger than $`1.2M_{}`$ don’t behave like the models. Hyades stars show a distinct “lithium dip” centered around $`1.4M_{}`$ Boesgaard & Tripicco (1986). Between $`1.2`$ and $`1.4M_{}`$ the lithium abundance is seen to decrease with increasing stellar mass. From $`1.4`$ to $`1.5M_{}`$ there is a sharp increase in lithium abundance. Balachandran (1995) finds the highest abundances for stars with masses above $`1.5M_{}`$ (the “blue side” of the dip) and for stars of $`1.2M_{}`$ (the “red side” of the dip); these stars have abundances about equal to those found in meteorites Grevesse & Anders (1989), abundances which are generally believed to be primordial.

Observations of young stellar clusters such as $`\alpha `$ Per and the Pleiades find little evidence of lithium depletion in more massive stars Soderblom et al. (1993). Indeed, there is a strong correlation between cluster age and Li abundance in the “lithium dip” between $`1.2`$ and $`1.5M_{}`$, indicating that the depletion takes place on the main sequence on a time scale of 100 Myrs or more Boesgaard (1991).

Observations of clusters having different metallicities $`[\mathrm{Fe}/\mathrm{H}]`$ show that the location of the dip is slightly metallicity dependent Balachandran (1995). Lower metallicity clusters have dips centered at lower masses. However, the high mass side of the dip depends only weakly on metallicity.

As we noted above, the lithium depletion appears to be due to mixing rather than gravitational settling or stellar mass loss. Two observations lead to this conclusion. First, surface depletions of lithium and beryllium are correlated, with the surface lithium diminishing more rapidly than surface beryllium Deliyannis et al. . Because the Be rich zone has twice the mass of the Li rich zone, mass loss would deplete all the Li before depleting the surface Be. Similarly, settling predicts that surface Li and Be would be depleted at roughly the same rate. Second, evolved stars with $`1.2<M<1.5`$ are also deficient in lithium Gilroy (1989); Balachandran (1995); settling of lithium would leave a non-depleted layer of lithium below the surface (in the absence of mixing), which would be dredged up to the surface when the convection zone deepens. The fact that the predicted increase in evolved stars is not seen indicates that the Li missing from main sequence stellar photospheres is destroyed rather than sequestered below the photosphere.

Above $`1.5M_{}`$ there appears to be little mixing of the upper layers of the star with deeper layers, since there is little Li depletion. However, the lithium in these stars does appear to be destroyed in deeper layers; evolved (giant) 1.6 solar mass stars are severely depleted Gilroy (1989). This depletion occurs inside the progenitors, below the photosphere Vauclair (1991), so that when the surface convection zone forms and then thickens, the surface lithium is mixed over a larger volume of depleted material, reducing the photospheric lithium abundance.

The surface mixing seen in the dip appears to be rotationally induced Charbonnel et al. (1992); Balachandran (1995), although this is far from certain.

As already noted, stars in $`\alpha `$ Per and the Pleiades clusters, with estimated ages of 50 and 70 million years, show little depletion in Lithium abundance for stars between $`1.2`$ and $`1.5`$ solar masses, while stars in the Hyades and older clusters show depletions that increase with cluster age. We interpret this to mean that the upper layers of the stars are mixed downward on time scales of $`100`$ million years.

The observational data summarized here leads us to propose the following simple model for the mass $`M_{mix}`$ of the surface mixing region in solar mass stars (see Figure 2). We first describe the lithium dip. The location, in mass, of the dip is described by $`M_{red}`$, $`M_{central}`$, and $`M_{blue}`$, giving the low, central, and high mass points of the dip. These three masses depend on the stellar $`[\mathrm{Fe}/\mathrm{H}]`$.

Above $`M_{blue}`$ we assume $`M_{mix}=M_{mh}`$; we use $`M_{mh}=3\times 10^3M_{}`$, the value obtained for the meridonal circulation model of Charbonnel et al (1992). This number is currently unconstrained by observation on the low side, but must be less than $`5\times 10^3M_{}`$ to prevent Li depletion.

For $`M_{central}M_{}<M_{blue}`$ we use a linear interpolation between $`M_{mh}`$ and $`M_{mc}3\times 10^2`$, where the latter value (for a star at the center of the dip) allows for complete Li depletion but only partial Be depletion. We also use a linear interpolation between $`M_{mc}`$ and $`M_{ml}`$ for $`M_{red}M_{}<M_{central}`$. Stars with $`M=M_{red}`$ appear to be only slightly depleted, so we take $`M_{ml}5\times 10^3M_{}`$.

Extrapolating from data in Balachandran (1995), we take

$`M_{red}/M_{}=`$ $`1.51+0.054[\mathrm{Fe}/\mathrm{H}]`$ (2)
$`M_{central}/M_{}=`$ $`1.42+0.25[\mathrm{Fe}/\mathrm{H}]`$ (3)
$`M_{blue}/M_{}=`$ $`1.31+0.31[\mathrm{Fe}/\mathrm{H}].`$ (4)
Note that the high mass or blue side of the dip is essentially independent of $`[\mathrm{Fe}/\mathrm{H}]`$. For stars with mass below $`M_{red}`$, we take the mass of the convection zone ($`M_{cvz}`$) plus the mass of the mixing zone $`M_{ml}`$ of a $`M_{red}`$ star on the low mass or red edge of the dip.

The mass of the convection zone at an age of $`10^8`$ years (the mass at $`10^9`$ is similar, see Figure 1) is found from the following empirical fit, derived from our stellar models:

$$M_{cvz}=a\left([\mathrm{Fe}/\mathrm{H}]\right)\times (M_3M_{})^{3.65},$$
(5)
where $`a\left([\mathrm{Fe}/\mathrm{H}]\right)=0.2[\mathrm{Fe}/\mathrm{H}]+0.49`$. The quantity $`M_3`$ is the (metallicity dependent) stellar mass at which the convection zone mass equals $`10^3M_{}`$,

$$M_3=a_3[\mathrm{Fe}/\mathrm{H}]+b_3,$$
(6)
where $`a_3=0.32`$ and $`b_3=1.45`$. Equation (5) is very accurate near $`M_3`$ and good to about $`20\%`$ at $`0.8M_{}`$.

Our model for the depth of the mixing region in a solar metallicity star is shown in Figure (2). Note that stars of different metallicity will have “lithium dips” (peaks in this figure) at slightly different masses than those shown in the figure.

Armed with this model, we can predict the average change in metallicity of a sample of stars subjected to iron accretion. Essentially we will look for behavior similar to that seen in the abundance of lithium. Polluted low mass stars of a fixed age will show an increase in metallicity with increasing mass (below the low mass end of the Li dip, around 1.0 solar masses). Intermediate mass polluted stars will show a decrease in \[Fe/H\] with increasing mass, up to $`M_{dip}`$, followed by a sharp increase at $`M_{blue}`$.

We want to stress that, while we are using the lithium data as a proxy to estimate the depth of the surface mixing layer of main sequence solar type stars, and while we anticipate that the surface abundances of lithium and iron will follow similar patterns, the reason for the abundance variations differ. The lithium abundances vary because lithimum is destroyed when it is mixed below the surface of the star. Stars in the mass range we consider neither destroy nor create iron. Rather, we assume that some iron is added to the surface of the star after it forms; it is then diluted by the mixing process down to the same depth that the lithium is mixed to. In other words the variations in lithium are produced by the destruction of lithium inside the star, while the variations in iron are produced by the addition of iron from outside the star.

## 4 OBSERVATIONAL CONSTRAINTS ON POLLUTION

In order to look for mass-dependent variations in \[Fe/H\], we must be able to estimate the masses of stars, and we must have reliable measurements of the metallicity. To find the mass of a star we need to know its luminosity, or alternately its absolute V magnitude $`M_v`$, its colors (we use $`BV`$), and its composition \[Fe/H\]. The HIPPARCOS catalogue gives parallaxes accurate to roughly one milliarcsecond for a selection of stars in the solar neighborhood, allowing us to find absolute magnitudes when we are given $`V`$. We use SIMBAD values for $`B`$ and $`V`$. We take metallicities from the (Cayrel de Strobel et al., 1997) catalogue. This lists spectroscopically determined values for \[Fe/H\] taken from the literature.

We take all the HD stars in the Cayrel de Strobel et al. catalogue for which there are HIPPARCOS parallaxes larger than 10 milliarcseconds (corresponding to distances less than 100pc). We also require that the error in the parallax be less than 10 percent. Since we are interested in main sequence stars that have masses less than about $`2M_{}`$, we eliminate all stars with $`M_v<1.0`$. We then fit the remaining stars to stellar tracks taken from a grid of models having $`0.5[\mathrm{Fe}/\mathrm{H}]0.5`$ in steps of $`0.05`$ dex, and masses in the range $`0.60M1.75M_{}`$ in steps of $`0.05M_{}`$. The models are from the Yale stellar evolution code.

A star is defined to be evolved if it has a convection zone mass more than ten times larger than the mass of the surface mixing layer at an age of $`10^8`$ years, while the absolute $`V`$ magnitude has not decreased by more than $`0.5`$. We refer to such stars as “Hertzsprung gap” stars. This occurs roughly when the star leaves the main sequence. It is a useful definition because it allows us another check on the pollution scenario; if pollution is occurring, the class of stars with convection zone masses ten times larger than their initial main sequence convection zone masses should show substantially lower \[Fe/H\] values than their parent populations.

If the absolute $`V`$ magnitude has decreased by more than $`0.5`$, the star is considered to be a giant. Because their surface gravities are so much lower than main sequence stars of the same mass, we hesitate to compare the metallicity trends of giants and dwarfs.

In the course of our analysis a problem with our stellar models became apparent. Stars with $`BV>1`$ were systematically found to be older than 20Gyrs. We believe that this is due to a failure of our atmospheric models at these low effective temperature. For this reason we do not include any stars redder than this in our analysis; this corresponds to stars with masses less than $`0.750.80M_{}`$. At the other end of the mass spectrum, we felt confident in extrapolating up to masses of $`1.80M_{}`$, $`0.05M_{}`$ beyond our grid.

Of the 650 stars in our sample, 466 were unevolved and had masses in the range $`0.8M1.8`$. Fifteen stars had convection zones 3-10 times larger than their main sequence value (at $`10^8`$yrs); we treated these as unevolved. Nineteen stars had convection zones more than ten times as massive as the main sequence values; these are our Hertsprung gap stars.

In Figure (3) we plot $`[\mathrm{Fe}/\mathrm{H}]`$ against stellar age, while Figure (4) shows $`[\mathrm{Fe}/\mathrm{H}]`$ plotted against stellar mass. The metallicity decreases with increasing age, and increases with increasing mass. We also note that the Hertzsprung gap stars tend to have lower metallicity than unevolved stars of the same age and mass. Figure (5) shows the histogram of \[Fe/H\] for our sample. It is well-fit by a gaussian with mean $`[\mathrm{Fe}/\mathrm{H}]=0.095`$.

The trend with stellar age is expected on the basis of chemical evolution; younger stars are constructed using gas that has been contaminated by material processed in the interiors of massive stars. The more recent the stellar birth, the higher the fraction of processed material. The trend is seen even more clearly when we bin the data by age, as seen in Figure (6). The slope of a least squares fit to the \[Fe/H\] versus log(age) curve is $`0.21`$.

The trend of increased metallicity with increased stellar mass might arise from the age-metallicity trend, combined with the short lifetimes of massive stars. Since the stars we consider are either unevolved or only slightly evolved, they must have ages less than or roughly equal to their main sequence lifetimes. Massive stars have short lifetimes, and hence low ages (see Figure 7). They will therefore be more metal rich, on average, than less massive stars.

In the next section we will show that, while part of the variation of metallicity with stellar mass is due to the age-metallicity relation, it is not the whole story.

## 5 MONTE CARLO MODELS OF POLLUTION

The data suggest a possible trend of surface \[Fe/H\] with stellar mass. In order to evaluate the significance of this trend, we have undertaken Monte Carlo experiments to calculate the predicted metallicity of a population of polluted stars.

We start by producing a population of unpolluted stars having the same characteristics as the Cayrel de Strobel-derived sample described above. This involves generating a sample of stars with mass roughly uniformly distributed between $`0.8`$ and $`1.8M_{}`$, with metallicities gaussian distributed about an age dependent mean. From our sample of stars and using a simple two parameter fit we find that the mean \[Fe/H\] is given by

$$[Fe/H]=\alpha \times \mathrm{log}(Age)+\beta ,$$
(7)
where $`\alpha =0.21\pm 0.02`$ and $`b=0.03`$. The width of the $`[\mathrm{Fe}/\mathrm{H}]`$ distribution is slightly age dependent, but we assume an age-independent width $`\sigma =0.18`$.

We constrain the stellar age to be less than the calculated lifetime of the star. As noted above we define the lifetime of a star to be the age at which the mass of the surface convection zone exceeds ten times the mass of the mixed surface layer at an age of 100 Myrs. Our stellar models have lifetimes (as defined above) which vary as the $`3.65`$ power of the stellar mass,

$$L(M,[\mathrm{Fe}/\mathrm{H}])=M_1([\mathrm{Fe}/\mathrm{H}])\times \left(\frac{M_{}}{M_{}}\right)^{3.65},$$
(8)
where

$$M_1([\mathrm{Fe}/\mathrm{H}])=11.5+4.5[\mathrm{Fe}/\mathrm{H}]Gyr$$
(9)
is the lifetime of a one solar mass star of metallicity \[Fe/H\].

Using this prescription we generate 466 stars, yielding the mass-age relationship depicted by the open circles in Figure 7; the agreement with the observed age-mass relationship is very good.

### 5.1 The Unpolluted \[Fe/H\]-Mass Correlation

We can use this Monte Carlo model to see if the \[Fe/H\]-mass correlation seen in our sample could be due simply to the known age-mass and age-\[Fe/H\] correlation. To repeat, the idea is as follows: massive stars are necessarily young. Young stars are metal rich, since they form out of gas that has been processed through many generations of high-mass (supernova producing) stars. Thus more massive stars should have higher metallicity than less massive stars.

This argument is partially borne out by the Monte Carlo model. We have adjusted the slope of the $`[\mathrm{Fe}/\mathrm{H}]`$ vs. logarithm of stellar age relation to obtain agreement with that seen in our sample ($`0.21`$ dex per log(Gyr)), then examined the metallicity as a function of stellar mass. The stellar metallicities do show a dependence on stellar mass. However, the slope of the $`[\mathrm{Fe}/\mathrm{H}]`$ vs mass relation for the observed stars is $`0.26\pm 0.03`$ dex per solar mass, while that of the unpolluted Monte Carlo model is $`0.18\pm 0.03`$, different by 3 standard deviations.

A second indication that there is a real mass dependence in the data comes from allowing for a third parameter in a least squares fit. Heretofore we have used a two parameter fit, an age-metallicity slope and an intercept. If we allow for a linear variation of \[Fe/H\] with mass, the slope of the fit is a third parameter. A fit of the form

$$[\mathrm{Fe}/\mathrm{H}]=a_1+a_2\mathrm{log}(age)+a_3\left(\frac{M}{M_{}}\right)$$
(10)
does give a significantly improved chi squared.

To see if this three parameter fit can really distinguish between a pure age dependence and a mixed age and mass dependence, we modified our Monte Carlo model to allow for a linear dependence of average metallicity on stellar mass. We then ran the resulting stars through our least squares routine. If we set the amount of added iron to zero but force an age dependence, the least squares fit finds the correct age-metallicity slope. It also correctly finds that there is no mass dependence, despite the apparent dependence in a metallicity versus mass plot. If we introduce a mass dependence into the model, the least squares fit correctly finds both the age and mass dependences. Finally, we produced samples of stars with a mass dependence but no age dependence. A plot of metallicity versus age clearly shows a trend; this is expected since we force the more massive stars, which are younger, to have higher metallicities. Despite this apparent age-metallicity trend, the linear least squares fit finds the correct mass-metallicity slope, and a zero slope for the age-metallicity relation.

These results suggest that the real data exhibits a mass trend independent of stellar age, though by itself the argument is not very convincing. More compelling evidence is obtained by comparing to more sophisticated models. The inadequacy of a simple age dependent metallicity is most clearly seen by examining the binned metallicity versus mass distribution, Figure (8). Another indication is the poor $`\chi ^2/dof`$ of the fit; the best unpolluted model, shown as the dotted line in the figure, gives $`\chi ^2=3.2`$ for the binned distribution. We proceed to examine more realistic polluted models, which give substantially better fits to the data.

### 5.2 Realistic Polluted Models

As noted above, pollution introduces a third free parameter, the amount of iron added to the star. In the last subsection we assumed that $`[\mathrm{Fe}/\mathrm{H}]`$ increases linearly with stellar mass. This is not what we expect from pollution, since the mass of the surface mixing layer is not linear, or even monotonic with stellar mass. In this section we assume a gaussian distribution of accreted mass, with a variance equal to half the mean. Stars allotted a negative amount of accreted mass by this process are assumed to have no added material. We assume the accreted mass is mixed over a surface layer of mass $`M_{mix}`$, as described in section 3.1. Our model for the mass of the surface mixing layer has many parameters, but we hold them fixed at the values given in that section. The single adjustable parameter is then the mean mass of accreted iron. We adjusted this mean to obtain a minimal $`\chi ^2`$ in the binned data, both metallicity-mass and metallicity-age. The resulting fit to the metallicity-mass data is shown as the solid line in Figure (8). We were able to find fits with reduced $`\chi ^2`$ equal to one.

The reduction in the $`\chi ^2`$ of the fit afforded by the more sophisticated pollution model is dramatic. The detailed variation of $`[\mathrm{Fe}/\mathrm{H}]`$ with stellar mass, including the steep increase starting at $`0.8M_{}`$, the dip around $`1.4M_{}`$, and the steep increase to $`1.5M_{}`$ seen in the binned data is reasonably similar to that seen in the model. We feel that the low $`\chi ^2`$ and the detailed variations of metallicity with mass offer strong support both for 1) a mass dependent mixing, as suggested by the Li dip; and 2) for accretion of $`0.5\pm 0.2M_{}`$ of iron in most stars in the solar neighborhood.

### 5.3 Hertzsprung-gap Stars

The metallicity-mass trend, combined with the jump in average \[Fe/H\] around $`1.5`$ solar masses strongly suggests that stars in our sample have accreted iron rich material. However, the properties of our sample are not well constrained, since the selection criteria of the various samples collated by Cayrel de Strobel are difficult to ascertain. Perhaps the apparent variations in average metallicity are due to some unknown selection effect. For example, one of the subsamples is that of Favata, Micela, & Sciortino (1997). A plot of \[Fe/H\] versus mass for this sample shows an increase toward large masses, starting at about $`1.15M_{}`$. However, this sample was selected to have $`BV>0.5`$; this eliminates low metallicity high mass stars, which are bluer than this limit. This is the origin of the increase in \[Fe/H\] at $`1.15M_{}`$ in this subsample.

We do not see an increase in \[Fe/H\] at $`1.15M_{}`$ in the large sample. However, perhaps the sharp jump we see at $`1.5M_{}`$ is due to a similar selection effect. The best way to eliminate this possibility would be to perform an unbiased survey of a large number of stars in the solar neighborhood. We encourage observers to undertake such a survey.

We do have available to us a small control sample. Stars that are slightly evolved, in the sense that their convection zones are 10 or more times larger than those on the main sequence, but which are not yet giants (having $`M_v`$ within $`0.5`$ of the main sequence value) provide a check on our interpretation. These stars sit in the Hertzsprung gap. Unlike giants, their surface gravities and fluxes are very similar to stars on the main sequence, so their metallicities can be compared directly to those of main sequence stars.

If pollution is responsible for the \[Fe/H\] variations over and above those due to age effects, than Hertzsprung-gap stars should have lower average metallicities, and they should show no trends with stellar mass.

The open triangles in Figure (4) represent the metallicities of stars in the gap. It is apparent that these stars have, on average, slightly lower metallicity than unevolved stars of similar mass. This is confirmed by an examination of the binned data, Figure (8), which also shows that there is no evidence for a change in average \[Fe/H\] in the gap stars at the blue edge of the lithium gap around 1.5 solar masses. This provides support for the notion that the variations in metallicity with mass seen in the unevolved stars is due to pollution.

The dashed line in Fig. (8) shows the predicted $`[\mathrm{Fe}/\mathrm{H}]`$ for the same model as the solid line, but where the mass of the accreted material has been reduced by a factor of ten. This is equivalent to increasing the mass of the surface mixing layer by the same factor. The predicted metallicity for this population is lower than that predicted for the polluted sample, and shows only a hint of a “lithium” dip. It gives an acceptable fit to the rather sparse observational data (reduced $`\chi ^2`$ less than one).

In plotting both the data and the fit we have not adjusted $`[\mathrm{Fe}/\mathrm{H}]`$ according to the star’s age. In Figure (9) we replot the data adjusting for the age-metallicity trend using the measured slope of the age-metallicity relation in Figure (6). This involves increasing the $`[\mathrm{Fe}/\mathrm{H}]`$ for old stars. This figure should be compared to plots of lithium abundance in clusters; we see evidence of a “lithium” dip in the iron data similar to that seen in lithium data. The dip is less distinct in the iron data, a fact we attribute to the range of metallicities in our sample (the metallicities in cluster stars show much smaller dispersions). The variation in metallicity produces a variation in the mass of the dip, so combining stars with different metallicity will tend to smear out the dip. This smearing is clearly seen in our Monte Carlo models (the solid line in Fig. 8).

We have also examine the metallicities of giants in our sample. They show a slight increase in \[Fe/H\] with mass, but this increase is consistent with the age-\[Fe/H\] correlation. The average \[Fe/H\] is actually slightly larger than either the unevolved stars or the Hertzsprung gap stars, although consistent within the rather large errors. However, given the very different environment for line formation in these stars, we do not feel it is appropriate to make a direct comparison between the two populations.

## 6 DISCUSSION

We made use of three free parameters in fitting the metallicity-age and metallicity-mass data; one for the slope of the metallicity-age relation, one (the mean of the added mass) for the correlation in the metallicity-mass relation, and an overall normalization. We want to stress that, in fitting the Hertzsprung gap star data, (both age and metallicity) we did not adjust any free parameters.

We also had a number of parameters in the model which were fixed by other observations; most were used to describe the location and depth of the lithium dip. The most critical for our purposes is the mass $`M_{mh}=3\times 10^3M_{}`$ of the surface mixing layer in stars blueward of the dip. Varying this mass will directly affect our estimates of the amount of iron accreted onto the star. For this reason our estimate of $`0.4M_{}`$ of accreted iron is uncertain. However, the need to accrete some material is not uncertain; reducing $`M_{mh}`$ will reduce the inferred mass of accreted material, but models with no pollution will still fail to fit the data.

The linear model for the mass of the surface mixing layer we employ is unlikely to be correct in detail; the “v” shaped bottom is particularly questionable. Future work should employ a more realistic model for $`M_{mix}`$, but for this first foray we did not feel that more realistic models were necessary.

If the trends we see are due to pollution, one might surmise that stars with larger \[Fe/H\] will form with more massive planetesimal disks. This might result in larger amounts of iron rich material falling on stars with larger intrinsic metallicities. This could be modeled by assuming a correlation between \[Fe/H\] and the amount of added material; once again we leave such experiments for later work.

### 6.1 Gas Migration and Ingestion of Gas Giant Planets

The gas giant planets in our solar system appear to have rock/ice cores with masses of order $`10M_{}`$, with the possible exception of Jupiter Guillot 1999a , which could lack such a core. These cores are believed to form at $`5`$ AU from the sun since the surface density of solid matter in the protoplanetary disk is believed to jump up beyond that distance due to the presence of ice. The discovery of Jupiter-mass planets in small (less than $`0.05`$ AU) orbits around solar type stars motivated a number of groups to suggest that gas giant planets migrate. If this migration is overly efficient, Jupiter-mass objects may be accreted onto the star. In the leading theory of planetary migration, that of tidal interactions between massive planets and the gas disk out of which they formed, the migration time is given by the viscous evolution time of the disk, prompting many authors to predict such accretion Lin (1997); Laughlin & Adams (1997). Others have suggested accretion of iron rich material originally between the planet and the star Gonzalez (1997).

There are two distinct scenarios. In the first, iron rich material is pushed onto the star by the gaseous disk. This material necessarily falls on the star before the gas disk disappears, when the star is between one and ten million years of age. In the second, the enriched material arrives via some other mechanism, possibly long after the gas disk disappears. The primary difference between the two, from the point of view of pollution, is the depth of the stellar convection zone; stars younger than ten million years and less massive than $`1.31.4M_{}`$ have very massive convection zones and consequently are difficult to pollute.

However, as pointed out by Laughlin and Adams (1997), stars more massive than $`1.31.4M_{}`$ will show evidence for pollution if they accrete iron rich material at ages $`10`$ Myrs. Figure (10) shows a Monte Carlo model where $`1.4M_{}`$ of iron is added to stars before they reach ten million years of age. We have adjusted the age-metallicity slope and the fraction of stars accreting material to achieve the best $`\chi ^2`$ fit; we can find acceptable fits (with $`\chi ^21`$) to the unevolved stellar data as illustrated by the solid line in the figure. Roughly $`1/4`$ of these massive stars must accrete Jupiter-analogues as they approach the main sequence.

This early accretion model does not do so well when compared to the evolved stars. The model predicts that the metallicities of the massive stars will lie below those of unevolved stars, as observed. However the model predicts that low mass stars will not show any reduction in $`[\mathrm{Fe}/\mathrm{H}]`$, nor is the drop predicted for the massive stars as large as that seen in our data. Using our best fit model (the solid line in Figure 10) and then increasing the mass of the surface mixing layer by a factor of ten gives the dashed line in the figure; for this model the reduced $`\chi ^2=2.5`$; the predicted metallicities are systematically higher than those observed.

There is also a theoretical difficulty with the early accretion model. The model is motivated by the possibility that gas in the protoplanetary accretion disk can push Jupiter-mass objects into the star. The difficulty is that the gas that does this is also likely to accrete onto the star. Since any excess iron sequestered in the doomed planet had to come from the gas in the disk, the gas in the disk must be iron poor. If this iron poor gas accretes at roughly the same time as the planet, the net change in the surface metallicity of the star will be zero.

The fact that the net change in metallicity of the star should be small, combined with the poor fit of the model to the evolved star data, suggest that accretion must occur after the accretion disk vanishes.

### 6.2 Ingestion of Gas Giants By Other Means

Lin (1997) speculates that Jupiter-mass planets might survive the gas disk only to fall onto the star later. If this occurs, and if these objects have heavy element abundances similar to gas giants in our solar system, then the star will gain $`1.8M_{}`$ of iron for each accreted object. We can find acceptable fits to the data in Fig. (8) if roughly $`30\%`$ of stars in the solar neighborhood have each accreted a single Jupiter-mass bodie after their convective envelopes thinned. These models naturally predict that the average metallicity of such stars will decrease by the observed amount when their convection zones deepen at the end of their main sequence lifetime.

However we have more information available to us than just the mean metallicity as a function of stellar mass. In particular, we know that the distribution of \[Fe/H\] is well fit by a gaussion, as illustrated in Figure (5). This fact allows us to rule out pollution dominated by late ingestion of Jupiter-mass planets having iron contents similar to that of the giant planets in our solar system. The idea is to use the sensitivity of the high mass stars to pollution. Adding $`1.8M_{}`$ of iron in discrete lumps to $`30\%`$ of the stars (to give the observed mean \[Fe/H\]) will produce a population with metallicities about three times that of the unpolluted population in these high mass stars. Using the Monte Carlo model that reproduces the variation of mean metallicity mentioned in the previous paragraph, we obtain the two-peaked histogram shown in Figure (11). This can be compared with the observed distribution for those stars above $`1.4M_{}`$ shown on the same figure. If these stars are in fact polluted, the pollution is distributed smoothly over most of the stars rather than being concentrated in a fraction of order $`30\%`$. We conclude that the metallicity trend we find is not produced by the late ingestion of Jupiter-mass objects.

### 6.3 Surface Metallicity Enhancement By Mass Loss or Dust Accretion

The solar wind is observed to be iron rich; in fact many elements with low first-ionization potentials are more abundant, while He, with its very high first-ionization potential, is under abundant Meyer (1993); Geiss (1998). If a substantial portion of the stellar convection zone were lost in a massive wind then the photospheric iron abundance would be less than the bulk abundance. If less massive stars suffer larger mass losses (relative to their more massive convection zones) this would produce a population in which less massive stars have lower metallicity than more massive stars, as we have found. When these low mass stars evolve, deep convection zones will form and mix metal rich material up to the surface; the high mass stars will also show a slight increase in iron abundance.

However, what we appear to see is that as massive stars evolve, their metallicities decrease; we don’t know what the low mass stars do when they evolve. There is a second difficulty with this scenario. The mass loss rates required of low mass stars are excessive; to lose a substantial fraction of the convection zone ($`0.05M_{}`$ for a $`0.8M_{}`$ star) the wind must carry away $`5\times 10^{11}M_{}/yr`$ for a billion years. The solar mass loss rate is $`3\times 10^{14}M_{}/yr`$. Recent limits on the mass loss rates in K dwarfs are much lower, $`10^{12}M_{}/yr`$ or less Lim & White (1996); van den Oord & Doyle (1997).

Of course the observation that lithium and beryllium depletions are correlated indicate that mass loss does not extend down to masses of order $`3\times 10^2M_{}`$ in stars of mass $`1.21.4M_{}`$. Taken together, these facts suggest that mass loss is not responsible for the metallicity variations we find.

Many young stars of roughly solar or larger mass show infrared excesses due to circumstellar dust Aumann et al. (1984); Habing et al. (1999); Decin et al. (2000). Much of this dust may end up on the star. Dust particles which absorb photons moving radially and then radiate them isotropically (in the dust rest frame) lose angular momentum in a process called the Poynting-Robertson effect. The amount of angular momentum lost ($`\delta \mathrm{}`$) per photon of frequency $`\nu `$ is $`\delta \mathrm{}/\mathrm{}h\nu /mc^2`$, where $`m`$ is the rest mass of the dust particle. Thus, for each accreted dust particle, a star must emit $`\mathrm{}/\delta \mathrm{}`$ photons, which means it must emit a total energy approximately equivalent to the rest mass energy of the accreted material. The total energy emitted during the period of dust accretion is larger than this by $`1/\tau `$, where $`\tau `$ is the optical depth to dust absorption. Thus, for a star radiating at approximately solar luminosity, the accretion of material will take a time

$$T_{dust}\frac{10^8years}{\tau }\left(\frac{M_{acc}}{2M_{}}\right)\left(\frac{L_{}}{L}\right)$$
(11)
For $`\tau 1`$, this number corresponds roughly to the observed lifetime of dust around young massive (F) stars Habing et al. (1999), but the spherically averaged $`\tau `$ is obviously considerably less than 1; values near $`10^4`$ are typical. Decin et al. (2000) find longer lived disks ($`1`$ Gyr) of similar optical depth, but the maximum amount of dust that such stars can accrete is still much less than an Earth mass. We conclude that dust accretion is unlikely to provide enough material to explain the apparent iron enhancements we see.

### 6.4 Accretion of the Interstellar Medium

Observations of the interstellar medium in the vicinity of the sun indicate that it has roughly solar abundances ($`+0.05`$ dex) of some elements, such as Zn, P, and S Howk et al. (1999). Iron is typically depleted, probably onto dust grains. The youngest stars in our sample have $`[\mathrm{Fe}/\mathrm{H}]0.2`$ (Figure 3), slightly larger than the ISM, but our results indicate that about $`0.050.1`$ of this is due to pollution. This is consistent with the notion that the bulk metallicity of young stars will equal that of the ISM out of which they form.

The implication is that older stars have metallicities lower than that of the ISM in which they currently reside. Accretion of the relatively metal rich ISM material could in principle produce metallicity signatures of the type we have found. However, it seems that the youngest stars in our sample may be metal rich compared to the present day ISM. In view of the rather large errors in the metallicity determinations we feel that this conclusion has to be considered as tentative.

### 6.5 $`\beta `$ Pictoris

The star $`\beta `$ Pictoris is surrounded by a dusty disk Smith & Terrile (1984). The dust is believed to have a life time shorter than the estimated age of the star (between $`3\times 10^7`$ and $`10^8`$ years) and so must be replenished, possibly by collisions between planetesimals Lecavelier Des Etangs et al. (1996). Transient red-shifted absorption features have been seen in spectra of $`\beta `$ Pic for the last decade Boggess et al. (1991); Lagrange et al. (1996). These events have been interpreted as infalling comets or asteroids. Events with very high velocities ($`>100`$km/s) are seen a few times per month Beust et al. (1991); Beust & Morbidelli (2000). Similar transients are seen around numerous other young stars de Winter et al. (1999); Grady et al. (1999).

The absorption features are believed to be produced by the cometary tails of the infalling bodies; Beust & Morbidelli argue that the infalling bodies must have radii larger than $`10`$ km in order to survive the evaporation long enough to evolve onto the inferred highly eccentric orbits. They estimate that $`1.8M_{}18M_{}`$ of material is removed from the planetesimal disk by this process if it has continued for $`10^8`$ years. They refer to evaporation of these bodies; however a small but substantial fraction of them will survive long enough to strike the star. This could provide enough material to give the signature we have seen.

### 6.6 Solar Pollution

The sun is metal rich for stars of its age and mass, as can be seen from Figs. (3) and (4). The implication is that the radiative interior of the sun might have a metallicity slightly lower than the photosphere. A number of authors Joss (1974); Christensen-Dalsgaard et al. (1979); Levy & Ruzmaikina (1994); Jeffery et al. (1997) suggested this possibility two decades ago as a way out of the solar neutrino problem; lowering the metallicity of the solar interior reduces the radiative gradient, which in turn would lower the inferred temperature in the core. Finally, a lower core temperature would predict that fewer neutrinos would be emitted compared to the standard solar model.

The metallicity deficit required to explain the lack of solar neutrinos is much larger than the $`0.017`$ dex increase resulting from the accretion of $`0.4M_{}`$ of iron. This avenue for solving the neutrino problem was dropped by solar workers, partly for this reason. More recently, Henney & Ulrich (1998) examined polluted models to see if accurate measurements of the the solar five minute oscillations (p-modes) could be used to detect differences between the metallicities of the surface convection zone and the radiative interior. They found that current data was unable to distinguish between unpolluted models and models with $`10M_{}`$ of accreted cometary material. The effect on the production of solar neutrinos was also negligible.

As an aside we note that the SOHO satellite has discovered a large number of sun grazing stars, a class of comets that plunge into the sun (see Raymond et al. (1998) and references therein). Roughly 200 have been observed by SOHO over its lifetime to date. Raymond et al. give an estimate of about $`10^3`$cm for the radius of the comet they observed. If most of the observed objects are of this size, and if we assume that the currently observed rate has been more or less steady over the lifetime of the sun, (one strikes the sun roughly every day as observed in 2000), then the mass of accreted material would be about one-one millionth of an Earth mass. Note however that most of the comets seen by SOHO appear on dynamical grounds to be the result of the break-up of a single object, so that the currently observed flux is likely to be a burst of transient activity.

## 7 CONCLUSIONS

We analyzed some 642 stars from the Cayrel de Strobel catalogue having spectroscopically determined metallicities, and well determined HIPPARCOS parallaxes. Using a large grid of stellar models, we determined the age and mass of the stars in our sample. We then examined the variation of metallicity with stellar age and mass.

We have found striking mass-dependent variations in photospheric iron abundances of main sequence solar mass ($`0.81.8M_{}`$) stars in the solar neighborhood. These variations mimic the variations seen in lithium abundances. With somewhat less confidence, because of the small sample size, we find that the iron abundances of Hertzsprung-gap stars are on average lower than those of main sequence stars, and that the metallicities of these slightly evolved stars have no mass dependence. These results are consistent with the accretion of an average of $`0.4M_{}`$ of iron onto the surface of main sequence stars. This strongly suggests that terrestrial-type material is common around solar type stars in the solar neighborhood.

We would like to thank Barth Netterfield for helpful conversations. Support for this work was provided by NSERC of Canada, and by NASA through Hubble Fellowship grant #HF-01120.01-99A from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. This research made use of the SIMBAD database, operated at CDS, Strasbourg, France.


text:
# Note on the Spin-Spin Relaxation Time in High-𝑇_c Cuprate Superconductors [1]

## Abstract

Recent experiments on the spin-spin relaxation time in high-$`T_\mathrm{c}`$ cuprate superconductors are discussed on the basis of the BCS-RPA theory. Keywords: spin-spin relaxation time, high-$`T_\mathrm{c}`$ cuprate superconductor

The temperature dependence of the spin-spin relaxation time recently measured in high-$`T_\mathrm{c}`$ cuprate superconductors  deviates from the expected behavior of a prevailing thoery. The observed relaxation time has strong temperature dependence for temperatures below the superconducting transition temparature, $`T<T_\mathrm{c}`$, while the theory predicts little temparature dependence. Such a temperature independence predicted by the theory was regarded as an evidence for the $`d`$-wave superconductivity and established by the measurements in La- and Y-systems. On the other hand, the measurements in Tl-, Hg- and Bi-systems  revealed the insufficiency of the theory.

In this Short Note we clarify that the temperature independence is not a direct consequence of the theory by searching much wider parameter space than in the original work and that the recent experiments can be understood within the framework of the theory.

The Gaussian component of the spin-spin relaxation time $`T_\mathrm{G}`$ is expressed as

$$\frac{T_G(T_c)}{T_G(T)}=\left[\frac{R(T)}{R(T_c)}\right]^{1/2},$$
(1)
using

$$R(T)=\underset{𝐪}{}F^4(𝐪)\chi ^2(𝐪)\left[\underset{𝐪}{}F^2(𝐪)\chi (𝐪)\right]^2.$$
(2)
Here $`\chi (𝐪)`$ is the static spin susceptibity and $`F(𝐪)`$ is the form factor.

The static spin susceptibity is given by the BCS-RPA theory as

$$\chi (𝐪)=\frac{\chi _0(𝐪)}{1U\chi _0(𝐪)},$$
(3)
with

$$\chi _0(𝐪)=\frac{1}{N}\underset{𝐪}{}\left[C_+(𝐩,𝐪)\frac{f(E_{𝐩+𝐪})f(E_𝐩)}{E_𝐩E_{𝐩+𝐪}}+C_{}(𝐩,𝐪)\frac{f(E_{𝐩+𝐪})f(E_𝐩)}{E_𝐩+E_{𝐩+𝐪}}\right],$$
(4)
where $`U`$ is the parameter for the on-site Coulomb interaction and $`f(E)`$ is the Fermi distribution function. The energy for the BCS quasiparticle is given by

$$E_𝐩=\sqrt{\epsilon _𝐩^2+\mathrm{\Delta }_𝐩^2},$$
(5)
with

$$\epsilon _𝐩=2t[\mathrm{cos}(p_xa)+\mathrm{cos}(p_ya)]4t^{}\mathrm{cos}(p_xa)\mathrm{cos}(p_ya)2t^{\prime \prime }[\mathrm{cos}(2p_xa)+\mathrm{cos}(2p_ya)]\mu ,$$
(6)
where $`t`$, $`t^{}`$, $`t^{\prime \prime }`$ is the hopping parameter and $`\mu `$ is the chemical potential and

$$\mathrm{\Delta }_𝐩=\frac{\mathrm{\Delta }_0(T)}{2}[\mathrm{cos}(p_xa)\mathrm{cos}(p_ya)],$$
(7)
where $`\mathrm{\Delta }_0(T)`$ is the gap function in the BCS theory. Here we consider a tight-binding model on the square lattice of the lattice constant $`a`$. The coherence factor is given by

$$C_\pm (𝐩,𝐪)=\frac{1}{2}\left[1\pm \frac{\epsilon _{𝐩+𝐪}\epsilon _𝐩+\mathrm{\Delta }_{𝐩+𝐪}\mathrm{\Delta }_𝐩}{E_{𝐩+𝐪}E_𝐩}\right].$$
(8)
The form factor is given by

$$F(𝐪)=A_{cc}+2B[\mathrm{cos}(q_xa)+\mathrm{cos}(q_ya)],$$
(9)
in the present two-dimensional tight-binding model.

We have done numerical calculations with above formulae. The parameters are fixed as $`t^{}/t=1/6`$, $`t^{\prime \prime }/t=1/5`$, $`\mu /t=1`$, $`T_\mathrm{c}/t=0.1`$ and $`U/t=3`$. The temperature dependence of the chemical potential $`\mu `$ has been neglected. The choice of $`T_\mathrm{c}/t=0.1`$ is in accordance with ref. 5 and $`U/t=3`$ is consistent with ref. 7. Our main results have been checked to be insensitive to the choice of the values of $`t^{}`$ and $`t^{\prime \prime }`$.

First, in Fig. 1(a), we have examined the temperature dependence of the static spin susceptibility. The susceptibility is strongly temperature dependent around $`𝐪=(0,0)`$ and weakly dependent around $`𝐪=(\pi /a,\pi /a)`$. The latter reflects the fact that the existence of the $`d`$-wave superconductivity does not conflict with the antiferromagnetic long-range correlation.

Since the form factor is fixed as $`A_{cc}=4B`$ in ref. 5, the temperature dependence of $`T_\mathrm{G}`$ is predominantly determined by that of $`\chi (𝐪)`$ around $`𝐪=(\pi /a,\pi /a)`$. Namely, the contribution from the region around $`𝐪=(0,0)`$ is masked by the form factor so that $`T_\mathrm{G}`$ has little temperature dependence. Such a chice is consistent with the data for La- and Y-systems. For example, $`A_{cc}/B=177/44`$ in YBa2Cu4O8.

On the other hand, the masking is absent in other cuparate superconductors. For example, $`A_{cc}/B=170/100`$ in Tl2Ba2CuO6 and the numerical result for $`T_\mathrm{G}`$ with this value in Fig. 1(c) has substantial temperature dependence. The 15% reduction at the lowest temparature is twice of the reduction in the theory of ref. 5 and consistent with the experiment.

It can be said that we have substantial temperature dependence if the masking is absent. However, the temperature dependence observed in Hg- and Bi-systems is much stronger than obtained in Fig. 1(c). As seen in Fig. 1(b) where $`A_{cc}/B=170/100`$, the shorter the coherence length of the superconductivity, the stronger the temperature dependence is. Although more elaborate treatment is needed for the short coherence-length (strong coupling) case, the qualitative feature can be captured by the present parameter extrapolation. In ref. 5 the parameter has been fixed as $`2\mathrm{\Delta }_0(0)/k_\mathrm{B}T\mathrm{c}=3.52`$ of the weak coupling value. Thus it can be said that we have a strong temperature dependence if the coherence length is short. Actually in the numerical result with parameters appropriate for Bi-system of Fig. 1(d), we have a strong temparature dependence as commented in ref. 4.

In summary, the temperature dependence of the spin-spin relaxation time observed by recent experiments can be understood within the BCS-RPA theory, if the coherence length of the superconductivity is short and the masking due to the form factor is absent. Such a conclusion is easily understood by taking account of the fact that the BCS susceptibility $`\chi _0(𝐪)`$ is suppressed by the Cooper pair (singlet) formation for $`|𝐪|<\xi ^1`$ where $`\xi `$ is the coherence length for the superconductivity (the size of the Cooper pair).

Since our calculation is based on the mean-field approximation, it can not directly fit experimantal data and will be modified if fluctuations are taken into account. The superconducting transition tempetatute $`T_\mathrm{c}`$ in the present theory corresponds to the temperature giving the maximum of $`1/T_\mathrm{G}`$ in experiments and is higher than experimentally observed $`T_\mathrm{c}`$. The steep decrease of $`1/T_\mathrm{G}`$ around $`T=T_\mathrm{c}`$ in the theory will become gentle due to fluctuations.

As noted in ref. 4, the on-set temperature for the decrease of $`1/T_\mathrm{G}`$ due to superconducting fluctuations in experiments is lower than the so-called pseudogap temperature observed by the spin-lattice relaxation time of the NMR and the ARPES experiments. The on-set temperature seems to correspond to the temperature giving the maximum of the Hall coefficient and this maximum temperature can be interpreted as the on-set temperature for superconducting fluctuations, since the collective contribution of superconducting fluctuations to the Hall coefficient has opposite sign to that of spin fluctuations dominant at higher temperatures.

The author is grateful to Kenji Ishida for valuable discussions on various NMR experiments.