text:
# THE MASSIVE DISK AROUND OH 231.8+4.2

## 1 INTRODUCTION

Solar-type stars typically lose over half their initial main sequence mass on the Asymptotic Giant Branch (AGB) before becoming white dwarfs (Habing 1996); mass loss is of central importance both in stellar evolution and the replenishment of interstellar matter. While many AGB stars exhibit approximately spherically symmetric mass loss which can be understood in terms of standard models of the interaction of stellar pulsations and radiation pressure on dust (Lamers & Cassinelli 1999, Willson 2000), others show markedly asymmetric outflows which are not well understood. Here, we report an investigation to learn more about the circumstellar material around OH 231.8+4.2, a star with a distinctive bipolar outflow.

OH 231.8+4.2 appears to be an evolved member of the open cluster M 46 (Jura & Morris 1985), and on the basis of the main sequence turnoff from this system, this mass-losing star probably had an initial main sequence mass of about 3 M. Both its cluster membership and “phase-lag” measurements in the reflection nebulosity (Kastner et al. 1992) indicate the distance to the star is approximately 1.3 kpc. Although the star is buried in an optically opaque circumstellar dust cloud, scattering in bipolar lobes perpendicular to the nominal disk indicates the presence of a luminous M9 star (Cohen & Frogel 1977). The star is variable with a period which might be as short as 648 days (Feast et al. 1983) or as long as 708 days (Kastner et al. 1992), and its luminosity ranges between 104 L and 2 $`\times `$ 104 L (Kastner et al. 1998). The luminosity and effective temperature can be explained if the star lies on the AGB.

At 60 $`\mu `$m, OH 231.8+4.2 is brighter than every other AGB star in the sky except IRC +10216 and the Egg Nebula (RAFGL 2688) (see Jura & Kleinmann 1989, Kleinmann et al. 1978). IRC+10216 is very bright at 60 $`\mu `$m because its distance is only $``$150 pc; it is not uniquely luminous. For OH 231.8+4.2, however, with an estimated distance of 1.3 kpc, then as found in the eDIRBE data from the COBE satellite, $`L_\nu `$(60 $`\mu `$m) can be as large as 3 $`\times `$ 1024 erg s-1 Hz-1 which is comparable to that from the Egg Nebula and, as shown in Figure 1, substantially larger than $`L_\nu `$(60 $`\mu `$m) displayed by most other mass-losing AGB stars.

The optical and near-infrared morphology of OH 231.8+4.2 indicates that there is a massive disk surrounding the central star (see, for example, Kastner al. 1992). Additionally, OH 231.8+4.2 displays a bipolar outflow at radio and optical wavelengths which is inferred to carry at least 0.3 M (Reipurth 1987, Morris et al. 1987, Sanchez-Contreras et al. 2000, Alcolea et al. 2001). The bipolar morphology exhibits discrete clumps (Alcolea et al. 2001).

Although not directly determined spectroscopically, in order to account for the bipolar morphology, a plausible assumption is that OH 231.8+4.2 is a binary and thus the outflow may be very flattened (Mastrodemos & Morris 1999, Soker & Rappaport 2000, Soker 2002). The nature of the putative companion is unknown. It is unlikely that the companion is a white dwarf since OH 231.8+4.2 does not exhibit any spectroscopic similarity to symbiotic systems.

The dynamical nature of the disk is not known. Is the matter expanding radially outwards from the star or is it gravitationally bound and orbiting the star, as appears to be the case for some other evolved red giants with circumstellar envelopes such as the Red Rectangle (Jura & Kahane 1999)? In order to learn more about OH 231.8+4.2, we have obtained high angular resolution mid-IR images of this system with the Keck I 10m telescope. Previously, Meixner et al. (1999) have obtained mid-IR images of this system, but with a 4m telescope. With the larger aperture of the Keck telescope, we can study more compact structures. Many other studies have focused on the large scale bipolar outflow; here we concentrate on the inner region of the system.

## 2 OBSERVATIONS

Our data were obtained on 2001 Feb 05 (UT) at the Keck I telescope using the Long Wavelength Spectrometer (LWS) which was built by a team led by B. Jones and is described on the Keck web page. The LWS is a 128 $`\times `$ 128 SiAs BIB array with a pixel scale at the Keck telescope of 0$`\stackrel{}{\mathrm{.}}`$08 and a total field of view of 10$`\stackrel{}{\mathrm{.}}`$2 $`\times `$ 10$`\stackrel{}{\mathrm{.}}`$2. We used the “chop-nod” mode of observing, and 2 different filters centered at 11.7 $`\mu `$m and 17.9 $`\mu `$m with widths of 1.0 $`\mu `$m and 2.0 $`\mu `$m, respectively. Following Chen & Jura (2001), we used Capella (= HR 1708) for flux and point-spread-function calibrations. For Capella the FWHM of the image was 0$`\stackrel{}{\mathrm{.}}`$47 and 0$`\stackrel{}{\mathrm{.}}`$49 at 11.7 $`\mu `$m and 17.9 $`\mu `$m, respectively.

The images in the 11.7 $`\mu `$m in the 17.9 $`\mu `$m filters are presented in Figures 2 and 3. The Keck data show an extension at position angle $``$22 that has a full length of $``$3″. The extended infrared emission which was reported by Meixner et al. (1999) shows the same orientation in the sky as the bipolar lobes detected at other wavelengths. There is also a central unresolved point-like source. We measure total fluxes of Fν(11.7 $`\mu `$m) = 27 Jy and Fν(17.9 $`\mu `$m) = 240 Jy. These fluxes are about 50% greater than the average values determined by IRAS of Fν(12 $`\mu `$m) = 19 Jy and Fν(25 $`\mu `$m) = 226 Jy. However, OH 231.8+4.2 is strongly variable at mid-infrared wavelengths. The DIRBE instrument on the COBE satellite measured mean fluxes of Fν(12 $`\mu `$m) = 72 Jy and Fν(25 $`\mu `$m) = 680 Jy. Thus our measured fluxes are within the range given by previous observations. If the pulsational period given by Kastner et al. (1992) has remained stable, then our data were obtained near minimum light. Both the COBE and IRAS beams are much larger than the angular resolution that we obtained, but the fluxes reported with these instruments may still be comparable to our results since we detect emission as far as 2$`\mathrm{}`$ from the central source. At this angular offset, we expect the grain temperature to have fallen to 170 K (see below) and thus there would probably be relatively little extended emission beyond this region at 11.7 $`\mu `$m and even 17.9 $`\mu `$m. We find that the unresolved central sources carry about 1/3 and 1/4 of the total flux at 11.7 $`\mu `$m and 17.9 $`\mu `$m, respectively, so that within the unresolved core, Fν(11.7 $`\mu `$m) = 9 Jy and Fν(17.9 $`\mu `$m) = 60 Jy.

It is notable that OH 231.8+4.2 does display a bright unresolved central source. Other stars with ratios of Fν(25 $`\mu `$m)/Fν(12 $`\mu `$m) $``$ 10 in the IRAS data, indicative of cold grains, show distinct shells and no central dust source. For example, HD 179821, with Fν(25 $`\mu `$m)/Fν(12 $`\mu `$m) = 21 shows a resolved shell with an inner diameter of about 3$`\stackrel{}{\mathrm{.}}`$3 (Jura & Werner 1999). In HD 179821, the weak 12 $`\mu `$m emission results from few grains being near the star. In the case of OH 238.1+4.2, there are a large number of grains near the star yet the material is relatively cold.

## 3 THE UNRESOLVED CENTRAL SOURCE

The infrared emission from OH 231.8+4.2 depends in part upon the composition and size of the grains. Here, we first consider the possibility that the grains are sufficiently small that the emissivity of the particles varies as $`\nu ^{+1}`$ and that the envelope is optically thin. In this case, the average temperature of the grains, $`T_{gr}`$ can be estimated implicitly from the formula:

$$\frac{F_\nu (11.7\mu m)}{F_\nu (17.9\mu m)}=\left(\frac{\nu _2}{\nu _1}\right)^4\frac{e^{h\nu _1/kT_{gr}}\mathrm{\hspace{0.17em}1}}{e^{h\nu _2/kT_{gr}}\mathrm{\hspace{0.17em}1}}$$
(1)
where $`\nu _2`$ and $`\nu _1`$ correspond to 11.7 $`\mu `$m and 17.9 $`\mu `$m, respectively. With $`\frac{F_\nu (11.7\mu m)}{F_\nu (17.9\mu m)}`$ = 0.15, then $`T_{gr}`$ $``$ 120 K. If the grains act like black bodies, then $`T_{gr}`$ $``$ 135 K is a better fit to the data.

If the material is opaque, then the observed flux, Fν, is given by the expression:

$$F_\nu =B_\nu (T_{gr})\mathrm{\Omega }_{source}$$
(2)
where $`B_\nu `$ denotes the Planck function and $`\mathrm{\Omega }_{source}`$ the solid angle of the source. Since the flux at 17.9 $`\mu `$m from the innermost 0$`\stackrel{}{\mathrm{.}}`$5 is 60 Jy, the minimum black body temperature of the grains which emit this radiation is 130 K.

In the simplest of optically thin models (see Sopka et al. 1985), if the grain opacity varies as $`\nu ^{+p}`$, then the temperature of the dust grains, $`T_{gr}`$, at a distance $`D`$ from the star can be computed from the expression:

$$T_{gr}=T_{}(\frac{R_{}}{2D})^{\frac{2}{4+p}}$$
(3)
Given that OH 231.8+4.2 has spectral type M9, we adopt $`T_{eff}`$ = 2500 K (van Belle et al. 1996). With an average value of $`L`$ = 1.5 $`\times `$ 104 L, then $`R_{}`$ = 4.6 $`\times `$ 1013 cm. Since the central source is unresolved, we expect that all the emission from the point source arises from a region within 0$`\stackrel{}{\mathrm{.}}`$25 of the star which corresponds to a distance of 4.9 $`\times `$ 1015 cm. At this distance, the predicted grain temperature from equation (3) for $`p`$ = 1 (a crude averaging of the silicate opacity given by David & Pegourie 1995) is 290 K, much greater than the value inferred from the data. Even if the particles are black bodies with $`p`$ = 0, then the grain temperature is predicted to equal 170 which is larger than the inferred value.

A better way to understand the data is to presume a nonspherical circumstellar envelope where much of the mass is largely confined to a disk. Knapp, Sandell & Robson (1993) showed that a flat opaque disk cannot explain the spectral energy distribution from OH 231.8+4.2 because in the infrared, such a disk should exhibit Fν varying as $`\nu ^{0.33}`$ while, in fact, between 12 $`\mu `$m and 60 $`\mu `$m, the IRAS data indicate that Fν varies as $`\nu ^{2.1}`$. Below, we consider a flared disk which, when viewed face-on, exhibits Fν varying as $`\nu ^{1.67}`$ (Chiang & Goldreich 1997). We suggest that emission from a flared disk at a substantial inclination can account for the data.

We can make an estimate of the minimum dust mass in the compact source, $`M_{dust}`$. If it were optically thin, then

$$M_{dust}=\frac{F_\nu D_{}^2}{B_\nu (T_{gr})\chi _\nu }$$
(4)
where $`D_{}`$ denotes the distance to the source. If, as is likely, the system is opaque, then this relationship only provides a lower bound to $`M_{dust}`$. With $`T_{gr}`$ = 130 K, Fν(17.9 $`\mu `$m) = 60 Jy and $`\chi _\nu `$ = 1000 cm2 g-1 (Ossenkopf, Henning & Mathis 1992), then $`M_{dust}`$ $``$ 7 $`\times `$ 1028 g.

## 4 A MODEL FOR THE DISK EMISSION

In order to explain the spectral energy distribution of the unresolved source, we present a model to explain the observations of OH 231.8+4.2 which is based on the calculations by Chiang & Goldreich (1997) for passive, opaque disks around pre-main sequence stars. Although we consider a post-main sequence system, the basic physical principles can be applied to both sorts of disks. In this model, the heating and cooling rates locally balance each other, and the task is to determine the temperature which achieves this effect. We ignore the illumination of the disk by light from the bipolar lobes.

If the disk has a local temperature, $`T_{disk}`$, then in a flat system which is optically thick in the vertical direction, for $`D`$ $`>>`$ $`R_{}`$,

$$T_{disk}\left(\frac{2}{3\pi }\right)^{1/4}\left(\frac{R_{}}{D}\right)^{3/4}T_{}$$
(5)
For a disk in vertical hydrostatic equilibrium which is locally isothermal, the density distribution, $`\rho `$, is given by the expression

$$\rho =\rho _0exp\left(\frac{z^2}{h^2}\right)$$
(6)
where

$$h=\left(\frac{2D^3k_BT_{disk}}{GM_{}\mu }\right)^{1/2}$$
(7)
and $`\mu `$ is the mean molecular weight. In this expression, we ignore the self-gravity of the disk.

The transition from flat to flared disk occurs at a critical distance, $`D_{crit}`$, which is established by the criterion that $`h`$ $``$ $`R_{}`$. From the above,

$$D_{crit}=\left(\frac{3\pi }{32}\right)^{1/9}\left(\frac{GM_{}\mu }{k_BT_{}R_{}}\right)^{4/9}R_{}$$
(8)
and at this location, the temperature, $`T_{crit}`$, is given by the expression:

$$T_{crit}=\left(\frac{256}{81\pi ^4}\right)^{1/12}\left(\frac{k_BT_{}R_{}}{GM_{}\mu }\right)^{1/3}T_{}$$
(9)
Using the stellar parameters for OH 231.8+4.2 described above, assuming that $`M_{}`$ = 1 M because the star has lost mass during its post-main sequence evolution, and assuming that the gas is primarily H2 and He so that $`\mu `$ = 3.9 $`\times `$ 10-24 g, then $`(2k_BT_{}R_{})/(GM_{}\mu )`$ = 0.061. Consequently, $`D_{crit}`$ = 3.0 $`R_{}`$ so that the flaring of the disk begins close to the star. At $`D`$ = $`D_{crit}`$, the disk temperature, $`T_{crit}`$, is 740 K.

In their models for flared disks, Chiang & Goldreich (1997) define the “grazing angle”, $`\alpha `$, as the angle between the surface of the disk and the line of sight to the star. Far from the star,

$$\alpha D\frac{d}{dD}\left(\frac{h}{D}\right)$$
(10)
In this approximation,

$$T_{disk}=\left(\frac{\alpha }{2}\right)^{1/4}\left(\frac{R_{}}{D}\right)^{1/2}T_{}$$
(11)
With these expressions, the disk temperature is given by the expression:

$$T_{disk}=\left(\frac{1}{7}\right)^{2/7}\left(\frac{R_{}}{D}\right)^{3/7}\left(\frac{2k_BT_{}R_{}}{GM_{}\mu }\right)^{1/7}T_{}$$
(12)
Chiang & Goldreich (1997) discuss models for pre-main sequence stars whose disks result from the collapse of interstellar clouds. Such disks might have very large amounts of angular momentum and extend to great distances from the star. Here, we picture a disk created in a binary system, and the initial angular momentum of the system limits the amount of angular momentum, $`J`$, that the disk can possess. Below, we argue that the disk has an angular momentum given by:

$$J=M_{disk}\sqrt{\frac{GM_{}D_{out}}{\pi }}$$
(13)
where $`D_{out}`$ is the outer boundary of the disk which we set equal to 3 $`v_{vis}t`$ where these parameters are defined below (see equation 18). If OH 231.8+4.2 had, itself, an initial mass of 3 M and a companion of 1 M in a circular orbit of radius between 3 and 5 AU, then the companion would have had an orbital angular momentum of $``$2 $`\times `$ 1053 g cm2 s-1. If half of this orbital angular momentum has been transferred to the disk so that $`J`$ $``$ 1053 g cm2 s-1, and if $`M_{disk}`$ $``$ 0.1 M, we find from equation (13) that $`D_{out}`$ = 6 $`\times `$ 1015 cm, consistent with the observational upper limit that $`D_{out}`$ $``$ 5 $`\times `$ 1015 cm. According to our models, where $`D_{out}`$ = 5 $`\times `$ 1015 cm, the inferred disk temperature from equation (12) is 130 K and $`h`$ = 3 $`\times `$ 1015 cm.

If the disk were nearly face-on-then we would detect emission from essentially every annular ring. However, if the angle between the normal to the disk and the line of sight, $`\theta `$, is greater than 45, then the inner, warm portions of the disk are occulted (Chiang & Goldreich 1999). According to Kastner et al. (1992) and Shure et al. (1995), $`\theta `$ $``$ 54, and therefore most of the disk is shadowed from our point of view. As a first approximation, we assume that we only observe emission from the outer portion of the flared disk which has radius, $`D_{out}`$, and a uniform temperature, $`T_{out}`$. This outer material is so “flared” into our line of sight that its emission dominates the unresolved central source that we detect. If $`\mathrm{\Omega }_{out}`$ denotes the solid angle subtended by the outer portion of the flared disk from our perspective, then

$$\mathrm{\Omega }_{out}\mathrm{\hspace{0.33em}4}\left(\frac{D_{out}}{D_{}}\right)^2sin\theta $$
(14)
We can use equation (2) to estimate the flux from the source.

The comparison of the model with the data is shown in Figure 4 where we present the results for both $`D_{out}`$ = 3 $`\times `$ 1015 cm (and $`T_{out}`$ = 130) and an alternative model with $`D_{out}`$ = 3 $`\times `$ 1015 cm (and $`T_{out}`$ = 160 K). The model with $`D_{out}`$ = 5 $`\times `$ 1015 cm agrees with our observations of the compact emission at 11.7 $`\mu `$m and 17.9 $`\mu `$m to better than a factor of 2. Our model does not reproduce the total flux from the system since our observations show that the bipolar lobes contribute the majority of the flux from OH 231.8+4.2 at 17.9 $`\mu `$m and 11.7 $`\mu `$m. The results from the model suggest that much of the millimeter and submillimeter continuum are produced by the compact source. This is possible but unproven. We do not know the millimeter-wavelength opacity of the particles around OH 231.8+4.2, and the disk may not be fully opaque at these wavelengths from the perspective of the Earth. In fact, the model significantly overestimates the flux at the longest observed wavelength of about 3 mm. The apparent agreement between the model and many of the integrated millimeter and submillimeter fluxes shown in Figure 4 may be a coincidence. Also, at 4.9 $`\mu `$m, the observed DIRBE flux is 16 Jy which agrees with ground-based measurements (Woodward et al. 1989) and, furthermore, this is the total flux expected from the photosphere. The map of Woodward et al. (1989) shows that the emission at 4.7 $`\mu `$m is extended and much of it is probably scattered radiation. The amount of disk thermal emission at wavelengths shortward of 11.7 $`\mu `$m is not easily measured and is therefore not included in Figure 4.

## 5 DISK DYNAMICS

The model proposed here is geometrically similar to that proposed by Cohen et al. (1985). What is uncertain is whether the disk material is expanding radially or whether it is orbiting; the expansion rate is not well measured. It has usually been assumed that the disk is expanding radially at the characteristic outflow speed of 20 km s-1. At this speed, the material reaches the inferred value of $`D_{out}`$ of 5 $`\times `$ 1015 cm in approximately 80 yr. If, however, the disk is in vertical hydrostatic equilibrium then the vertical sound crossing time must be short compared to the expansion time. In the outermost region, the disk temperature may be $``$100 K and therefore possess a sound speed near 1 km s-1. Since the vertical displacement approaches 5 $`\times `$ 1015 cm, the time to achieve vertical hydrostatic equilibrium and thus the implied disk lifetime is at least $``$2 $`\times `$ 103 years. It is easier to understand the presence of a flared disk if the material is orbiting.

Pringle (1991) has presented a model for how viscous disks around binary stars may expand as angular momentum is transfered from the binary stars into the disk, and here we follow his model. In the standard treatment of accretion disks, matter flows inwards while angular momentum flows outwards. Pringle (1991) argues that for a circumbinary disk, the binary system can provide enough torque to the disk so that no matter flows inwards, but the outward flow of angular momentum still occurs. In this model, after the disk is formed, its mass remains constant. The torque operates as the binary creates a tide in the circumbinary disk (Lin & Papaloizou 1979).

We assume that the circumstellar matter is injected into the system as a narrow ring of mass, $`M_{disk}`$, and an inner radius, $`D_{in}`$. The ring then expands as angular momentum is transferred from inner to outer regions by viscosity, $`\nu `$, which can be written as

$$\nu =\frac{2\alpha c_S^2}{3\mathrm{\Omega }}$$
(15)
Here, $`\alpha `$ is the usual dimensionless parameter and is taken to be $``$ 1, $`c_S`$ is the speed of sound or ($`\gamma k_BT_{disk}/\mu `$)1/2 where $`\gamma `$ is the ratio of specific heats, and $`\mathrm{\Omega }`$ is the local angular velocity of the disk, $`(GM_{})^{1/2}D^{3/2}`$.

Following Pringle (1991), we make the approximation that the gas temperature falls as $`D^{1/2}`$ so that

$$T_{disk}=T_{in}\left(\frac{D_{in}}{D}\right)^{1/2}$$
(16)
where $`T_{in}`$ is the temperature at the inner boundary of the disk, $`D_{in}`$. This radial variation of the temperature is somewhat different from that described in equation (12), but it should not lead to any essential change in the physics of the outflow. It does allow for a simple analytic solution to the asymptotic disk evolution. Following Pringle, with this radial variation of the temperature, we define $`v_{vis}`$ \[denoted $`k`$ by Pringle\] as a characteristic speed of the outflow induced by the viscosity, such that $`v_{vis}`$ = $`\nu /D`$, and is given by the expression:

$$v_{vis}=\frac{2\gamma k_BT\alpha }{3\mu }\left(\frac{D}{GM_{}}\right)^{1/2}$$
(17)
Because $`T`$ scales as $`D^{1/2}`$, $`v_{vis}`$ is independent of $`D`$.

The inner radius of the disk might be $``$1.7 times greater than the orbital separation of the binary pair (Artymowicz & Lubow 1994). Although uncertain, we adopt $`D_{in}`$ = 3 $`R_{}`$ = 1.3 $`\times `$ 1014 cm where $`T_{disk}`$ = 740 K (see above). The value of $`\alpha `$ is not well known for the material around OH 231.8+4.2. Saslaw (1978) has shown that grains around a variable star move in non-circular orbits, and these non-circular motions might produce a substantial viscosity. We assume that $`\alpha `$ = 1 implying that $`v_{vis}`$ = 0.29 km s-1.

After the initial conditions damp out, Pringle’s model for the outer portion of the disk has a surface density given by the expression:

$$\mathrm{\Sigma }\frac{M_{disk}}{(12\pi ^3D^3v_{vis}t)^{1/2}}exp\left(\frac{D}{3v_{vis}t}\right)$$
(18)
where $`M_{disk}`$ is the mass of the disk. This expression for the surface density of the disk appears to be independent of the torque being exerted on the disk because the transfer of angular momentum from the binary to the disk is effected by employing the inner boundary condition of the solution for the disk evolution that no matter flows inwards even though matter and angular momentum flow outwards. That is, at all times, the mass of the disk is constant and:

$$_0^{\mathrm{}}2\pi D\mathrm{\Sigma }(D)𝑑D=M_{disk}$$
(19)
The disk angular momentum, $`J`$, is given by:

$$J(t)=M_{disk}\sqrt{\frac{3GM_{}v_{vis}t}{\pi }}$$
(20)
From equation (18), we write that $`D_{out}`$ $``$ 3 $`v_{vis}t`$.

According to Pringle (1991), the characteristic time for the torque to operate on the disk, $`t_{torque}`$, at the innermost radius, is given by the expression:

$$t_{torque}=\frac{4D_{in}^2}{3\nu _{in}}=\frac{2\mu \sqrt{GM_{}D_{in}}}{\alpha \gamma k_BT_{in}}$$
(21)
Thus, for the disk around OH 231.8+4.2, this characteristic time for the torque to operate is $``$ 200 years. The outflow time is $`D_{out}/v_{vis}`$ or about 5000 yr. Therefore, the time for the torque to be effected is considerably less than the inferred outflow time. Alcolea et al. (2001) derive a minimum age of the bipolar outflow from OH 231.8+4.2 of 770 years. Thus the time noted here for the flared disk to expand to its current estimated size is comparable to the minimum lifetime of the system.

## 6 DISCUSSION

We have suggested that the circumstellar disk around OH 231.8+4.2 is gravitationally bound and orbiting the system rather than simply undergoing a radial expansion. Underlying this picture is that OH 231.8+4.2 is an AGB star with a hitherto-unseen companion. Ignace, Cassinelli & Bjorkman (1996) have suggested that mass loss from a red giant with a sufficiently large rate of rotation – perhaps induced by a companion (Livio & Soker 1988) – could produce a disk material within a few stellar radii of the mass-losing star. Whether this disk is bound or unbound – that is, whether it is orbiting or expanding – is not easily calculated because the acceleration law of the material near the star is not known. In fact, the acceleration may be time-variable over several pulsational cycles (Winters et al. 2000). Although it is possible that an orbiting ring could be formed near a rotating AGB star, this is unproven. Once such a ring is formed, as discussed above, Pringle (1991) has proposed that the material could expand by the viscous dissipation of angular momentum contained within the orbiting central binary.

Although there are a large number of studies of the molecular gas around OH 231.8+4.2, there is no direct kinematic demonstration from observations of the circumstellar gas that the disk around OH 231.8+4.2 is orbiting. Most of the gas emission is dominated by material in the bipolar flow which subtends a solid angle of $``$ 10-8 ster. The flared disk probably only subtends a solid angle of $``$ 10-11 ster. Therefore, the disk emission is mostly overwhelmed by the extended emission from the bipolar flow, and in order to detect the disk, interferometry is required.

One molecule which has been extensively studied at high angular resolution is OH. A remarkable aspect of OH 231.8+4.2 is that the centroid of the OH maser emission is offset about 1″ to the Southwest from the SiO and H2O masers and the millimeter continuum (Gomez & Rodriguez 2001, Zijlstra et al. 2001, Sanchez-Contreras et al. 1998). One possibility is that since OH is produced from the photodissociation of H2O (Glassgold 1996), a one-sided ultraviolet flux might produce a one-sided concentration of OH, as might explain the CO profile toward HD 188037 (Jura et al. 1997). Since M46 with its concentration of ultraviolet-emitting A-type stars also lies to the Southwest of OH 231.8+4.6, we suggest that the displacement of the OH masers from the core of OH 231.8+4.6 may be a consequence of its position in the outer, Northeast portion of M46.

Alcolea et al. (2001) report interferometric maps of the CO emission with resolution as high as 1$`\stackrel{}{\mathrm{.}}`$5 $`\times `$ 0$`\stackrel{}{\mathrm{.}}`$7. They distinguish between the extended bipolar flow and a central core. Even the central core displays both the north-south outflow plus a central equatorial disk or torus. The innermost disk or torus is not well mapped in their data, and it is not possible to determine whether it is orbiting or expanding.

Sanchez Contreras et al. (2000) report interferometric observations of HCO+, SO, H13CN and SiO around OH 231.8+4.2. The SO is particularly sensitive to high density regions, but the angular resolution of the maps 4$`\stackrel{}{\mathrm{.}}`$2 $`\times `$ 2$`\stackrel{}{\mathrm{.}}`$2 is insufficient to distinguish an orbiting disk of radius 5 $`\times `$ 1015 cm. Sanchez Contreras et al. (2000) report an expanding disk or ring of radius 2 $`\times `$ 1016 cm at speed of 6-7 km s-1. The relationship between this observed large disk of 2 $`\times `$ 1016 cm radius and the orbiting disk of 5 $`\times `$ 1015 cm radius that we propose to exist is unclear. Interferometric observations of the circumstellar H2O masers report two velocity peaks separated by 0$`\stackrel{}{\mathrm{.}}`$1 in the North-South direction (Gomez & Rodriguez 2001). This orientation suggests that this emission may be excited in the bipolar flow rather than the disk.

On a scale of 0$`\stackrel{}{\mathrm{.}}`$005, Desmurs et al. (2001) report SiO maser emission which is in the East-West direction and with kinematics which are consistent with being infalling material at a speed of $``$10 km s-1 with some additional rotation at a speed, $`V_{rot}`$, of $``$6 km s-1. The inferred rotation is measured over a spatial scale of approximately the diameter of the star. Unless the system is nearly face-on, the apparent rotation is too small to be orbital motion; it must reflect rotation of the star. The total angular momentum of an AGB star is about $`\frac{2}{9}M_{}R_{}V_{rot}`$ (Soker 1998); OH 231.8+4.2 may currently possess $``$ 1052 g cm2 s-1 of angular momentum. Above, we suggest that the disk may possess $``$ 1053 g cm2 s-1 of angular momentum. Thus, it is imaginable that there might be enough angular momentum in the system to account for the proposed disk.

We conclude that the available interferometric data are consistent with the hypothesis of an orbiting disk with an outer radius of 5 $`\times `$ 1015 cm, but they do not lend any support to the model. The argument for the existence of an orbiting rather than expanding disk is indirect. While the predicted infrared emission from a flared, orbiting disk agrees with our data, we cannot rule out other models for the compact, mid-infrared emission.

The mass of the disk around OH 231.8+4.2 is uncertain. Above we derived a minimum mass of dust of 1029 g which would correspond to a total of about 0.01 M if the gas to dust ratio is $``$160. The maximum possible dust mass is about 1.0 M since this would require that most of the material lost by OH 231.8+4.2 was trapped in the disk. However, a disk of this mass could not be both orbiting and as large as 5 1015 cm because it would possess too much angular momentum. It is possible that the disk around OH 231.8+4.2 has about 0.1 M. In this case, the disk mass around OH 231.8+4.2 is larger than that for disks around other evolved red giants except, perhaps, the Red Rectangle (see Jura & Kahane 1999). With time, the circumbinary disk around OH 231.8+4.2 must evolve. If, in fact, the material is orbiting, then macroscopic solids and even planets might eventually form.

## 7 CONCLUSIONS

We have obtained 11.7 $`\mu `$m and 17.9 $`\mu `$m images of OH 231.8+4.2; we identify an unresolved central source with Fν(17.9 $`\mu `$m) = 60 Jy. We propose that this emission results from a flared disk with an outer temperature of 130 K and an outer disk radius of 5 $`\times `$ 1015 cm. One possible model that is consistent with the data is that this flared disk is orbiting rather than simply expanding.

This work has been partly supported by NASA. We thank J. Alcolea, E. Chiang, M. Morris and N. Soker for very helpful comments.

FIGURE CAPTIONS

Fig. 1. A histogram of $`L_\nu `$(60 $`\mu `$m) for the “very-dusty” AGB stars within $``$1 kpc of the Sun listed by Jura & Kleinmann (1989). The distances are taken from that paper and we mostly used IRAS data for the 60 $`\mu `$m fluxes. For both the Egg Nebula and OH 231.8+4.2, we take $`L_\nu `$(60 $`\mu `$m) from the time-averaged COBE data.
Fig. 2. The 11.7 $`\mu `$m image of OH 231.8+4.2. North is up and East is to the left. The contour levels (Jy arcsec-2) are shown in the color bar.
Fig. 3. The 17.9 $`\mu `$m image of OH 231.8+4.2. North is up and East is to the left. The contour levels (Jy arcsec-2) are shown in the color bar. The very faint emission to the South of the image is an artifact.
Fig. 4. The spectral energy distribution of OH 231.8+4.2 compared with the model described in the text for an opaque, flared disk. The color of the squares denotes their sources: blue (IRAS), red (COBE), cyan (IRAM, Sanchez Contreras et al. 1998, 2000), green (JCMT, Knapp et al. 1993) and black (Kuiper Airborne Observatory, Harvey et al. 1991). The crosses show our results reported in this paper; black for the fluxes from the unresolved point source and red for the fluxes from the entire source. The solid black line shows the model described in the text which is comprised of the emission from a black body at 130 K from a cylinder of radius 5 $`\times `$ 1015 cm at a tilt angle of 54. The dashed line shows an alternative model with a cylinder at 160 K and a radius of 3 $`\times `$ 1015 cm at the same tilt angle of 54.


format: text

text:
# Effects of stoichiometry, purity, etching and distilling on resistance of MgB2 pellets and wire segments

## 1 Introduction

Within weeks of the announcement of the discovery of superconductivity in MgB2 by Akamitsu and co-workers , it was established that high purity, very low residual resistivity samples of MgB2 could be synthesized by exposing boron powder or filaments to Mg vapor at temperatures at or near $`950^{}C`$ for as little as two hours . Samples with residual resistivity ratio \[$`RRR=R(300K)/R(42K)`$\] values in excess of 20 and residual resistivities as low as 0.4 $`\mu \mathrm{\Omega }cm`$ were synthesized by this method. Such a low resistivity in an intermetallic compound with a superconducting critical temperature, $`T_c`$, near $`40K`$ was of profound physical, as well as engineering, interest. The implications of this high $`RRR`$ and low $`\rho _0`$ ranged from large magneto-resistances (in accordance with Kohler’s rule) to questions of how a material with such an apparently large electron-phonon coupling could have such a small resistivity. On the applied side, a normal state resistivity of 0.4 $`\mu \mathrm{\Omega }cm`$ for temperatures just above $`T_c`$ means that MgB2 wires would be able to handle a quench with much greater ease than, for example, Nb3Sn wires which have a $`\rho _0`$ that is over an order of magnitude larger for $`T20K`$ .

Unfortunately other techniques of synthesizing MgB2 have not yet been able to achieve such high $`RRR`$ or low $`\rho _0`$ values . In some cases the authors of these papers have concluded that the resistivity of their samples must be the intrinsic resistivity and that higher $`RRR`$ values or lower residual resistivity values must somehow be extrinsic. In order to address these concerns and in order to shed some light on how low resistivity samples can be grown we have studied the effects of boron purity and magnesium stoichiometry on sintered pellet samples. In addition we have studied the effects of filament diameter and post synthesis etching and distilling on MgB2 wire segments. Based on these measurements we conclude that the purity of the boron used to make the MgB2 is a dominant factor in determining the ultimate, low temperature, normal state resistivity of the sample, and that $`RRR`$ values as high as 20 and residual resistivities as low as 0.4 $`\mu \mathrm{\Omega }cm`$ are intrinsic materials properties of high purity MgB2.

## 2 Sample synthesis

Samples of MgB2 for this study were made in the form of sintered pellets as well as wire segments. The sintered pellets were made by sealing stoichiometric amounts of Mg and B into Ta tubes and placing these tubes (sealed in quartz) into furnaces heated to $`950^{}C`$ for 3 hours and then quenched to room temperature . For the initial studies of boron purity stoichiometric MgB2 was synthesized and the quality of the boron was varied. For the studies of magnesium stoichiometry nominal stoichiometries that ranged from Mg0.811B2 to Mg1.211B2 were used and samples were synthesized with isotopically 99.95% enriched 11B. For all synthesizes lump Mg of 99.9% purity was used.

Wire segments of MgB2 were made by sealing boron filaments purchased from Textron  or Goodfellow  into a Ta tube with excess Mg, using a ratio of approximately Mg3B. After reacting the filaments for 68 hours at $`950^{}C`$ the Ta tubes were quenched to room temperature and the wire segments were removed from the Ta reaction vessel. Given that there can be some excess Mg on the surface of the wire segments, some of the wire segments were etched in a solution of 5% HCl and ethyl alcohol for times up to 5 minutes. This treatment removes the surface Mg and leads to the surface of the wire segments having the same appearance as the surface of the stoichiometric sintered pellets: a slightly golden / grey color. Another method was used to remove any potential surface Mg from the MgB2 wire segments: high temperature distillation of the Mg. In order to achieve this a wire sample was place into a quartz tube that was continuously pumped by a turbomolecular pump to a pressure less than $`10^5Torr`$. The evacuated tube was then heated to $`600^{}C`$ for 12 hours. This temperature and time were chosen in part because attempts at higher temperature distillation lead to a decomposition of the MgB2 in the wire itself.

A.C. electrical resistance measurements were made using Quantum Design MPMS and PPMS units. Platinum wires for standard four-probe configuration with connected to the sample with Epotek H20E silver epoxy. LR 400 and LR 700 A.C. resistance bridges were used to measure the resistivity when the MPMS units were used to provide the temperature environment. Powder X-ray diffraction measurements were made using a Cu $`K_\alpha `$ radiation in a Scintag diffractometer and a Si standard was used for all runs. The Si lines have been removed from the X-ray diffraction data, leading to apparent gaps in the powder X-ray spectra.

## 3 Effects of Boron Purity

Table 1: Boron form and purity (As provided by the seller).    Purity Form Source Main Impurities    90%    Amorphous  (325 mesh)  Alfa Aesar    Mg  5%   95%    Amorphous  ($`<`$ 5 mesh)  Alfa Aesar    Mg  1%   98%    Crystalline  (325 mesh)  Alfa Aesar    C  0.55%   99.95%    Isotopically pure 11B  Crystalline (325 mesh)  Eagle-Picher    Si  0.04%   99.99%    Amorphous  (325 mesh)  Alfa Aesar    Metallic Impurities  0.005%

Table 1 presents the source and purity information available for each of the starting boron powders used. Figure 1 presents powder X-ray diffraction spectra for three samples with varying nominal boron purities: 90% purity, 99.99% purity, and the 99.95% purity, isotopically pure 11B. By comparing the two upper panels to the bottom panel it can be seen that the strongest MgB2 lines are present in all three samples. The upper panel, the data taken on the sample made from boron with only a nominal 90% purity, also has a weak Mg and MgO powder lines present. This is not inconsistent with the fact that the primary impurity in the 90% boron is associated with Mg.

Figure 2 displays the normalized resistance, $`R(T)/R(300K)`$, of MgB2 pellets that were made using five different types of boron powder. Each curve is the average of three resistivity curves taken on different pieces broken off of each pellet. Figure 2 demonstrates that $`RRR`$ values can range from as low as 4 to as high as 20 depending upon what source of boron is used. Among the natural boron samples examined there is a steady increase in $`RRR`$ as the purity of the source boron is improved. The MgB2 synthesized from the isotopically pure boron appears to have the best $`RRR`$, although it’s nominal purity is somewhat less than that of the 99.99% pure natural boron, but those skilled in the art will realize that claims of purity from different companies can vary dramatically. In addition, it is very likely that the isotopically pure boron was prepared in a somewhat different manner from the other boron powders used (very likely using a boron fluoride or boric acid or any of its complexes as an intermediate phase in order to achieve isotopic separation). The primary point that figure 2 establishes is that the purity of the boron used can make a profound difference on the normal state transport properties.

In Fig. 3 the same resistance data is plotted, but instead of simply normalizing the data at room temperature the data is normalized to the temperature derivative at room temperature. This is done to see if the resistance curves differ only by a temperature independent residual resistivity term: i.e. this normalization is based upon the assumption that the slope of the temperature dependent resistivity at room temperature should be dominated by phonon scattering and therefore be the same for each of these samples. As can be seen this seems to be the case, at least to the first order. By using higher purity boron we are able to diminish the additive, residual resistance by a factor of approximately five.

The insets to figs. 2 and 3 also indicate that there is a monotonic improvement in $`T_c`$ as the boron purity (or $`RRR`$ value) is increased. $`T_c`$ values vary from just below $`38K`$ to just above $`39K`$ depending upon which boron is used.

Based upon these results we choose the isotopically pure 11B for the further study of the effects of Mg stoichiometry on MgB2 pellet samples. But before we proceed to the next section it is worth noting that one of the difficulties associated with the samples made by other research groups may well be due to the use of boron with less than the highest purity. In addition, to our knowledge very few other groups have been using the Eagle-Picher isotopically pure boron in the samples for electrical transport measurements.

## 4 Effects of Magnesium Stoichiometry

In order to study the effect of magnesium stoichiometry on the transport properties of Mg11B2 a series of Mgx11B2 ($`0.8x1.2`$) samples were synthesized. Figure 4 presents powder X-ray diffraction patterns for the extreme members of the series (top and bottom panels) as well as for the stoichiometric Mg11B2 (middle panel). In all cases the lines associated with the Mg11B2 phase are present. For the Mg0.811B2 sample there is a weak line seen at $`2\theta =35.8^{}`$ that is associated with MgB4 (marked with a +). This is consistent with the fact that there was insufficient Mg present to form single phase Mg11B2. For the Mg1.211B2 sample there are very strong diffraction lines associated with Mg (marked with \*). This too is consistent with the stoichiometry of the sample: Mg11B2 is the most Mg-rich member of the binary phase diagram, therefore any excess Mg will show up as unreacted Mg. The X-ray diffraction pattern for the stoichiometric Mg11B2 shows much smaller peaks associated with a small amount of both MgB4 and Mg phases. This pattern is different from the one shown in Fig. 2 in that this sample was reacted for three hours whereas the sample used in Fig. 2 was reacted for four hours. Given that all of the samples used for the Mg-stoichiometry study were reacted for 3 hours it is appropriate to show this powder diffraction set along with the other members of the series. It should be noted that there is continuous change in the nature of the second phases in the samples. For Mg deficient samples there is only MgB4 as a second phase. For the stoichiometric Mg11B2 samples there are either no second phases or very small amounts of both MgB4 and Mg (depending upon reaction times), and for the excess Mg samples there is no MgB4, but clear evidence of excess Mg.

Figure 5 presents normalized resistance data for eight different Mgx11B2 pellets. In each case the curve plotted is the average of three or more samples cut from the same pellet. There is far less variation between the different pellets in this case than there was for the case of boron purity (Fig. 2). Figure 6 plots the $`RRR`$ values for each of the individual samples (shown as the smaller symbols) as well as the $`RRR`$ of the average curve. As can be seen the $`RRR`$ values increase slowly from $`14`$ for Mg0.811B2 to $`18`$ for Mg11B2. This is followed by a clear increase in $`RRR`$ values for excess Mg, with Mg1.211B2 having an $`RRR`$ value of $`24`$. The important point to note is that even for the most Mg deficient sample the lowest measured $`RRR`$ value is significantly greater than 10. At no point in this series we find samples with $`RRR`$ values of 3, 6, or 10, even when a clear MgB4 second phase is present. For samples ranging from Mg0.911B2 to Mg1.111B2 (dotted box in figure 6) the average $`RRR`$ values cluster around $`RRR=18\pm 3`$. These data indicate that for sintered pellets $`RRR`$ values of 18 can be associated with stoichiometric Mg11B2.

Whereas the effects of excess Mg are relatively minor in these samples (given their low intrinsic resistivities) these effects can be clearly seen. In addition to the increase in the $`RRR`$ value there is a change in the form of the temperature dependence of the resistivity. This can be best seen in Fig. 7 in which the resistance data have been normalized its room temperature slope. The data for all $`x`$ values less than 1.0 are similar and collapse into a single curve. On the other hand the resistance data for the $`x=1.1`$ and $`x=1.2`$ are qualitatively different. They start out with somewhat higher normalized resistance data than the stoichiometric sample and then below $`100K`$ cross below the stoichiometric sample (Fig. 7 inset). This is very likely due to the increasing effects of having Mg in parallel (and series) with the MgB2 grains. As can be seen in Fig. 7 this effect becomes larger as the amount of excess Mg is increased. This deviation from the MgB2 resistance curve may actually serve as a diagnostic for the detection of excess Mg.

## 5 MgB2 wire segments: Effects of diameter, etching and distilling

MgB2 can also be synthesized in the form of wire segments . The starting material is boron filament that has a small tungsten-boride core. Upon exposure to Mg vapor the boron filament is transformed into MgB2 wire segments, a process that is accompanied by an expansion of the wire diameter. Table 2 presents data on the initial and final diameters of the boron fiber and MgB2 wire segments used for this study. There is a clear expansion associated with the transformation of the boron into MgB2. The average increase in the diameter associated with this reaction is $`1.4`$ times. There is some uncertainty associated with this number due to the fact that once the MgB2 is formed the wire segments have a somewhat irregular surface as well as variation of the diameter along the length of a segment. The tungsten boride core does not manifest a noticeable change in diameter during this process. The size of the tungsten boride core of the wire segments is listed on the far side of Table 2. For the rest of this paper the wire segments will be identified by the diameter of the initial boron filament used to create them.

Table 2: Main properties of MgB2 wires and resistivity at $`300K`$ and $`RRR`$ for as-grown, etched and distilled MgB2 filaments with different diameters.(Note: Samples by the diameter of initial boron filaments used.)       Approximate Diameter \[$`\pm 10\mu m`$\]     Initial (Boron) Final (MgB2)      Approximate  Expansion     WBx  Approximate  Diameter \[$`\mu m`$\]      100 140  1.4 15     140 190  1.5 15     190 290  1.5 15     300 370  1.3 20  Sample Wires $`\rho (\mathrm{𝟑𝟎𝟎}𝐊)[\mu 𝛀\mathrm{𝐜𝐦}]`$ $`\mathrm{𝐑𝐑𝐑}`$   100 13.0 18  AS 140 10.8 45  GROWN 190 8.3 25   300 14.4 36   100 17.5 39  ETCHED 140 15.8 30   190 9.1 45   300 16.3 28  DISTILLED 300 9.3 35

Given that the wire segments are synthesized in a Mg rich vapor (the nominal stoichiometry is Mg3B), there is concern that small amounts of excess Mg vapor condense onto the surface of the MgB2 wire segments during the cooling process. This could lead to contributions to the temperature dependent resistivity from metallic Mg. In addition there is the possibility that the tungsten boride core may act as a low resistance resistor in parallel with the MgB2. Measurements of the transport properties of MgB2 wire segments of varying diameters will allow us to examine, and ultimately discount, both of these concerns.

If there were to be a significant contribution from metallic Mg on the surface of the wire segments, and if it is to be assumed that the metallic Mg has a lower resistivity than the intrinsic resistivity of the MgB2 (an assumption that is supported experimentally by our earlier data on the Mg1.211B2), then the effect of this excess Mg would scale with the surface area to volume ratio of the sample: i.e. there would be a substantially larger effect seen for the smaller diameter wires than for the larger diameter wires. In a similar manner the potential effect of the tungsten boride core would scale with the square of the ratio of the tungsten boride diameter to the MgB2 wire diameter. Given that the tungsten boride diameter remains between 15 and $`20\mu m`$ over the whole series and that the MgB2 diameter increased from $`140\mu m`$ to $`370\mu m`$, the potential effect that the tungsten boride would have would be largest in the smaller diameter wires and smaller in the larger diameter wires.

The temperature dependencies of the normalized resistivity of MgB2 wire segments are shown in Fig. 8a. All four diameter wires have similar temperature dependencies, but manifest somewhat different $`RRR`$ values. The inset to Fig. 8a shows the low temperature behavior near $`T_c`$. From this plot it becomes clear that there is no apparent effect of metallic Mg or tungsten boride on the resistivity. The $`RRR`$ values for the 100, 140, 190 and $`300\mu m`$ wire segments are 18, 45, 41 and 25 respectively. The highest $`RRR`$ value is for $`140\mu m`$ and the lowest value is for $`100\mu m`$. To first order, there appears to be little or no correlation between the wire diameter and the $`RRR`$ values. If any trend is to be extracted it is that $`RRR`$ values are generally higher for the larger diameter wires, a trend that contradicts the assumption that metallic Mg or tungsten boride are affecting the resistivity measurements.

In order to further examine the possible effects of metallic Mg on the transport properties of the wires we etched the as-grown wire segments in a solution of 5% HCl in ethyl alcohol for 5 minutes. This lead to an apparent removal of any Mg coating on the wire surface. The temperature dependence of the normalized resistance of these etched samples are shown in Fig. 8b. In this case the values of $`RRR`$ are increased. If the excess Mg in wires were acting as a parallel resistance in the sense of shorting MgB2, we would expect that with its removal we would find smaller $`RRR`$ values, not larger.

In the distillation process, the MgB2 as-grown wire with $`300\mu m`$ boron fiber initial diameter, which has $`RRR=36`$ was heated to $`600^{}C`$ for 12 hours and submitted to a continuously pumped high vacuum for removal of any excess Mg. This process practically did not alter the value of $`RRR`$ 35. This is strong evidence that the high values of $`RRR`$ which we obtained for MgB2 wires are intrinsic and not an influence of excess Mg.

Table 2 presents our estimates of the room temperature resistivity for each wire sample. The average room temperature resistivity for these nine samples is roughly 13 $`\mu \mathrm{\Omega }cm`$. It should be noted that the acid etching does increase the measured resistivity of the wires. We believe that this is due to two effects: i) a reduction of cross sectional area that has not been accounted for and ii) the creation of small cracks in the wire. In both cases the will change the effective geometry of the sample.

A final point about all of the wire samples is worth noting: all of the measured $`RRR`$ values are comparable to or better than those found for the stoichiometric 11B pellet samples. These boron fiber are made using a boron fluoride intermediate step and are reported to be 99.999% pure. This again points out that the purity (and probably the purification process) of the boron may be a critical variable.

## 6 Conclusion

In summary, through the synthesis of various pellets of MgB2 with different types of boron we found values of $`RRR`$ from 4 to 20, which covers almost all values found in literature. To obtain high values of $`RRR`$, high purity reagents are necessary. With the isotopically pure boron we obtained the highest $`RRR`$ 20 for the stoichiometric compound. We also investigated Mgx11B2 samples with 0.8 $`<x<`$ 1.2. These have shown that from the most Mg deficient samples we observe inclusions of the MgB4 phase, and no evidence of Mg. For samples with excess Mg we do not observe any MgB4. For the range Mg0.811B2 up to Mg1.211B2 we found average values of $`RRR`$ between 14 and 24. For smaller variations in stoichiometry ($`x\pm 0.1`$) $`RRR=18\pm 3`$. In addition our study of MgB2 wires as function of diameter is consistent with pellet results and inconsistent with either Mg or tungsten boride core acting as resistor in parallel with MgB2 filaments. All of our data are point to conclusion that high $`RRR`$ ($`20`$) and low $`\rho _0`$ ($`0.4\mu \mathrm{\Omega }cm`$) are intrinsic materials properties associated with high purity MgB2.

Acknowledgements

Ames Laboratory is operated for the US Department of Energy by Iowa State University under Contract No. W-7405-Eng-82. This work was supported by the Director for Energy Research, Office of Basic Energy Sciences and the National Science Foundation under grant No. DMR-9624778. The authors would like to thank D. K. Finnemore, M. A. Avila and N. Anderson for helpful assistance and many fruitful discussions.


text:
# The Anderson Transition in Two-Dimensional Systems with Spin-Orbit Coupling

## Abstract

We report a numerical investigation of the Anderson transition in two-dimensional systems with spin-orbit coupling. An accurate estimate of the critical exponent $`\nu `$ for the divergence of the localization length in this universality class has to our knowledge not been reported in the literature. Here we analyse the SU(2) model. We find that for this model corrections to scaling due to irrelevant scaling variables may be neglected permitting an accurate estimate of the exponent $`\nu =2.73\pm 0.02`$.

According to the scaling theory of localisation Abrahams et al. (1979), if interactions are neglected, all states are localised in two dimensions (2D). Two exceptions predicted by the scaling theory to this oft recited mantra are the extended states which occur at the center of a Landau level in the quantum Hall effect (QHE) Huckestein (1995), and the Anderson transition which occurs in zero magnetic field if there is a significant spin-orbit interaction Hikami et al. (1980); Ando (1989).

A surprising aspect of the QHE is that the estimate $`\nu =2.35\pm .03`$ Huckestein and Kramer (1990); Huckestein (1995) for non-interacting electrons of the critical exponent $`\nu `$, which describes the divergence of the localisation length $`\xi `$ at the transition, is in close agreement with the measured value. This is despite the fact that interactions are clearly relevant since the dynamical exponent $`z`$ is predicted to be two for non-interacting electrons while the experimental value is unity.

Critical phenomena are determined by the symmetry of the Hamiltonian and the dimensionality of the system. The important symmetries for the Anderson transition exhibited by non-interacting electrons are time reversal symmetry and spin rotation symmetry. There are three universality classes: orthogonal, unitary and symplectic. Systems with time reversal symmetry but where spin rotation symmetry is broken by the spin-orbit interaction belong to the symplectic class.

In this paper we estimate using numerical simulation and finite size scaling the exponent $`\nu `$ for the Anderson transition in the symplectic universality class in 2D. In contrast to the QHE, where extended states occur only in vanishingly small energy region, the metallic phase extends over a finite energy interval. This system is, therefore, a good candidate for the study of 2D quantum phase transitions. Early work suggested that metallic phase in this model is destroyed when interactions between electrons are taken into account Altshuler and Aronov (1983); Castellani et al. (1984), while more recent work suggests that this is not so Castellani et al. (1998).

Recently, there have been numerous reports of the observation of a zero magnetic field metal-insulator transition in 2D together with measurements of the critical exponents which characterize this transition Abrahams et al. (2001). Whether or not those experiments indicate the existence of a true metallic phase at zero temperature remains in dispute. If there is indeed a transition, the physics driving it and whether there is a relation with the transition we study here is not yet clear. However, just as was the case for the QHE, we believe that an accurate estimate for the metal-insulator transition in zero field in non-interacting 2D systems may prove useful.

There has been only very limited success in estimating the critical exponent with field theoretic methods Hikami (1992); Brezin and Hikami (1997). In Table 1 we tabulate the estimates of the exponent reported in previous numerical studies Ando (1989); Fastenrath et al. (1991); Evangelou (1995); Schweitzer and Zharekeshev (1997); Merkt et al. (1998); Minakuchi (1998); Yakubo and Ono (1998). There is considerable variation between these estimates. The estimates reported in Fastenrath et al. (1991); Evangelou (1995); Merkt et al. (1998); Minakuchi (1998); Yakubo and Ono (1998) seem to be consistent with a true value of the exponent in the range $`[2.6,2.9]`$. However, the estimate of Schweitzer and Zharekeshev (1997) is somewhat below this and that of Ando (1989) is in contradiction with the estimates of Fastenrath et al. (1991); Minakuchi (1998).

We use the transfer matrix method Pichard and Sarma (1981); MacKinnon and Kramer (1983) to estimate the localization length $`\lambda `$ of electrons on an $`L\times L_z`$ quasi-1D strip and then extract the critical exponent from a finite size scaling analysis of the dependence of $`\lambda `$ on $`L`$ and disorder. Two important factors limiting the accuracy of the estimate of the exponent obtained in this way are the accuracy of the data for $`\lambda `$ and the maximum width $`L`$ for which data are available. The standard error in the estimate of $`\lambda `$ decreases as $`\sqrt{\lambda /L_z}`$, while for a fixed $`L_z`$ the CPU time needed increases as $`L^3`$. At the critical point $`\lambda `$ increases linearly with $`L`$ so that the CPU time needed to estimate $`\lambda `$ to a given accuracy increases as $`L^4`$. This means that it is somewhat easier to improve the accuracy of the numerical data than to increase the size of the systems simulated. However, as the accuracy of the raw data improves, corrections to scaling due to irrelevant scaling variables become more important. While such corrections can be taken into account Slevin and Ohtsuki (1999), the number of fitting parameters is increased and correspondingly the uncertainty in the estimate of the exponent is increased. It is therefore advantageous to choose a model for which such corrections are negligible even when the raw data are of high accuracy. In this paper we report results for an SU(2) model for which this condition is satisfied.

The Hamiltonian for the SU(2) model describes non-interacting electrons on a 2D square lattice with nearest neighbour hopping

$$H=\underset{i,\sigma }{}ϵ_ic_{i,\sigma }^{}c_{i,\sigma }V\underset{<i,j>,\sigma ,\sigma ^{}}{}R(i,j)_{\sigma \sigma ^{}}c_{i,\sigma }^{}c_{j,\sigma ^{}}$$
(1)
Here $`c_{i,\sigma }^{}`$ ($`c_{i,\sigma }`$) denotes the creation (annihilation) operator of an electron at the site $`i`$ with spin $`\sigma `$ and $`ϵ_i`$ denotes the random potential at site $`i`$. We assume a box distribution with each $`ϵ_i`$ uniformly and independently distributed on the interval $`[W/2,W/2]`$. The width $`W`$ of the distribution measures the strength of the randomness. The constant $`V`$ is taken to be the unit of energy, $`V=1`$.

The spin-orbit coupling appears in the hopping matrix $`R(i,j)`$ between each pair of nearest neighbours on the lattice. These matrices belong to the group SU(2) of 2 $`\times `$ 2 unitary matrices with determinant one. The hopping matrices are parameterised as follows

$`R(i,j)=(\begin{array}{cc}e^{i\alpha _{i,j}}\mathrm{cos}\beta _{i,j}& e^{i\gamma _{i,j}}\mathrm{sin}\beta _{i,j}\\ e^{i\gamma _{i,j}}\mathrm{sin}\beta _{i,j}& e^{i\alpha _{i,j}}\mathrm{cos}\beta _{i,j}\end{array})`$ (4)
This matrix describes a rotation of the electron spin in three dimensional space. (The Euler angles of this rotation are related to, but not equal to, the angles $`\alpha `$, $`\beta `$ and $`\gamma `$.) In the SU(2) model the distribution of these angles is chosen so that the $`R(i,j)`$ are uniformly distributed with respect to the group invariant measure (Haar measure) on SU(2). This corresponds to $`\alpha `$ and $`\gamma `$ uniformly distributed in the range $`[0,2\pi )`$, and $`\beta `$ distributed according to the probability density,

$$P(\beta )\mathrm{d}\beta =\{\begin{array}{cc}\mathrm{sin}(2\beta )\mathrm{d}\beta \hfill & 0\beta \frac{\pi }{2}\hfill \\ 0\hfill & \mathrm{otherwise}.\hfill \end{array}$$
(5)
Hopping matrices on different bonds of the lattice are statistically independent. Periodic boundary conditions are imposed in the transverse direction. The necessary calculations are carried out using quaternion arithmetic Ando (1989) which halves the required number of multiplications compared with an implementation using complex arithmetic.

Some of the physics of the SU(2) model can be understood by comparing it with the Ando model which has been adopted in Ando (1989); Fastenrath et al. (1991); Schweitzer and Zharekeshev (1997); Yakubo and Ono (1998). In the Ando model as the electron propagates through the material its spin precesses at a rate and about an axis which depend on the electrons wave number. When scattered by the random potential, the rate and the axis about which the electron’s spin rotates changes. This leads to a diffusive motion of the spin with an associated spin relaxation length. Quantum interference between time reversed electron trajectories longer than this length produces the weak anti-localization effect Bergmann (1982). Motivated by the conjecture that the spin relaxation length might be an important irrelevant length scale, we adopted the SU(2) model where the uniform distribution of the hopping matrices on SU(2) ensures that the spin relaxation length is the shortest possible. Doing so we do indeed find that corrections due to irrelevant scaling variables can be neglected.

To determine the critical exponent $`\nu `$, critical disorder $`W_c`$ and other critical properties of the transition we fit dependence of $`\mathrm{\Lambda }=\lambda /L`$ on the system size $`L`$ and the disorder $`W`$, or when $`W=0`$ on the Fermi energy $`E_F`$, to a one parameter scaling law of the form

$$\mathrm{ln}\mathrm{\Lambda }=F\left(\psi \right).$$
(6)
Here $`\psi `$ is the relevant scaling variable. We expand the scaling function as a power series

$$F\left(x\right)=\mathrm{ln}\mathrm{\Lambda }_c+x+a_2x^2+\mathrm{}+a_{n_0}x^{n_0}$$
(7)
terminating the expansions at order $`n_0`$. To allow for non-linearity of the scaling variable, the scaling variable is approximated by an expansion in terms of the dimensionless disorder $`w=(W_cW)/W_c`$, where $`W_c`$ is the critical disorder separating the insulating and metallic phases. (If $`W=0`$ we set $`w=(E_FE_c)/E_c`$.) The growth of the relevant scaling variable with system size is described by the critical exponent $`\nu `$

$$\psi =L^{1/\nu }\left(\psi _1w+\psi _2w^2+\mathrm{}+\psi _{n_\psi }w^{n_\psi }\right),$$
(8)
where we terminate the expansion at order $`n_\psi `$. This same exponent describes the divergence of the localization (correlation) length

$$\xi =\xi _\pm \left|\psi _1w+\psi _2w^2+\mathrm{}+\psi _{n_\psi }w^{n_\psi }\right|^\nu .$$
(9)
where we terminate the expansion at order $`n_\psi `$. The absolute scales $`\xi _\pm `$ of the localization length on either side of the transition are not determined in this analysis, so we set them both to unity for simplicity. The linear coefficient in the expansions of $`F`$ is set to unity, as shown, to eliminate some redundancy in the definition of the fitting parameters. The total number of parameters is $`N_p=n_0+n_\psi +2`$.

The best fit is determined by minimizing the $`\chi ^2`$ statistic. The quality of the fit is assessed with the goodness of fit probability $`Q`$. Confidence intervals for the fitted parameters are estimated using a Monte Carlo method Press et al. (1992). This involves using the model and the best estimates of the fitting parameters to generate an ideal data set. From this data set a large ensemble of synthetic data sets is generated by adding random errors, with a variance equal to that of the error of the corresponding data point, to the ideal data set. Fitting of the ensemble of synthetic data sets produces a distribution for the critical parameters from which confidence intervals and the goodness of fit are estimated. This procedure is standard and systematic but does not take into account any unknown systematic effects that might only be discernible for very much larger systems. Of course, this caveat applies to almost any numerical estimate of a critical exponent.

Before turning to the estimate of the critical exponent we sketch the phase diagram in Figure 1. The figure is based on data for systems with sizes ranging from $`L=8`$ to $`L=32`$ for a number of Fermi energies between $`E_F=0`$ and $`E_F=3`$ for which the critical disorder $`W_c`$ was estimated, and on data for $`W=0`$ for which the critical energy $`E_c`$ was estimated. In the absence of a random potential i.e. when $`W=0`$, the Hamiltonian may have chiral symmetry and, in addition to the transition at finite $`E_c`$, a critical state may also be present at the band center Brouwer et al. (2000). Whether this possibility is realised depends on the boundary conditions and whether the number of sites is even or odd. Since chiral symmetry is broken by a random potential it does not affect our estimation of $`\nu `$ below.

To estimate the critical exponent accurately more extensive simulations were carried out for a single energy $`E_F=1`$. The numerical data are presented in Figures 2 and 3. Data with an accuracy of $`0.1\%`$ are available for system sizes $`L=8,16`$ and $`32`$, with accuracy $`0.2\%`$ for $`L=64`$, and $`0.4\%`$ for $`L=96`$. This required $`L_z`$ of the order of $`10^7`$ to $`10^8`$ depending on the size, the disorder and the accuracy. When fitting the data the intervals of $`W`$ and $`\mathrm{ln}\mathrm{\Lambda }`$ to consider must be decided. The exact choice is not particularly important provided all data are in the critical regime. A larger interval of $`\mathrm{ln}\mathrm{\Lambda }`$ requires a higher order of expansion in Eq. 7, while a larger interval of $`W`$ requires a higher order expansion in Eq. 8. The results of the finite size analysis are presented in Tables 2, 3 and 4. A number of fits of the numerical data are possible but, as can be seen by referring to Table 2, all yield consistent results. The estimates of the critical parameters are also stable against restriction of the system sizes under consideration, see Table 3, and against a narrowing of the range of disorder, see Table 4. The lines shown in Figure 2 and 3 correspond to a fit with $`n_0=3`$ and $`n_\psi =2`$. To demonstrate single parameter scaling graphically we re-plot the data as a function of $`L/\xi `$ in Figure 3; metallic and insulating branches are clearly visible. These are described by related scaling functions $`F_+`$ and $`F_{}`$ where

$$F_\pm \left(x\right)=\mathrm{ln}\mathrm{\Lambda }_c\pm x^{1/\nu }+a_2x^{2/\nu }+\mathrm{}+\left(\pm 1\right)^{n_0}a_{n_0}x^{n_0/\nu }$$
(10)
In summary, we have studied the Anderson transition in the 2D symplectic universality class and estimated the critical exponent $`\nu `$. We find $`\nu =2.73\pm .02`$, where the error is a $`95\%`$ confidence interval. Our result is consistent with the estimates reported in Fastenrath et al. (1991); Evangelou (1995); Merkt et al. (1998); Minakuchi (1998); Yakubo and Ono (1998) but not with those of Ando (1989); Schweitzer and Zharekeshev (1997). Analyses based on energy level statistics, such as Schweitzer and Zharekeshev (1997), also seem to have a tendency to produce estimates which are lower than those of the transfer matrix method for the 3D orthogonal universality class Ndawana et al. (2002). On the other hand, in our opinion, the error bar for the exponent claimed in Ando (1989) is too optimistic.


text:
# Black hole entropy without brick walls

## Acknowledgements

The author thanks X. J.Wang, Z. Q. Bai and H. G. Luo for their zealous helps during the research. The author also thank Profs. Y. Z. Zhang, R. G. Cai and C. G. Huang for their comments on this research. This work is supported by the Post Doctor Foundation of China and K. C. Wong Education Foundation, Hong Kong.


text:
# An Approximate Light Cone Method to Investigate Meson Structure

## 1 Introduction

It is well-established that the light-front formulation of QCD make calculations of observables in high-energy phenomena like DIS considerably easier while providing a more natural interpretation of the same . Transverse lattice formulation of QCD attempts to exploit the benefits  from using the LF gauge and the well-known advantages of the more conventional Euclidean and Hamiltonian lattice formulations.

Transverse lattice techniques were employed quite successfully to pure glue calculations . Encouraged by these results, the femtoworm approximation, which we describe briefly in the next section, was conceived to include the fermions on the transverse lattice, as well.

In this paper, we briefly describe the femtoworm approximation in section 2, referring the reader to relevant references for the most part. Some typical results and conclusions are included in section 3.

## 2 The Hamiltonian of the Femtoworm Approximation

Using light cone gauge requires that the longitudinal light cone space, $`x^+`$ and light cone time, $`x^{}`$ variables be redefined in terms of the equal time longitudinal space, $`x^0`$ and time, $`x^3`$ coordinates as $`x^\pm =x^0\pm x^3`$. The transverse directions are the same as in equal time definitions but discretized. This enables the use of a gauge invariant regulator in the theory while simultaneously retaining boost invariance in the longitudinal direction. In addition, the femtoworm approximation is based on the following considerations:
$``$ We work in the large $`N_c`$ limit.
$``$ A light cone Tamm-Dancoff truncation is applied so that the Fock state expansion of the meson states include only 2-particle and 3-particle Fock components.
$``$ One link approximation: The maximum number of link fields that can separate a $`q`$ and an $`\overline{q}`$ is one.
$``$ Interactions between 2-particle to 3-particle states (and vice versa) can take place only through two interaction terms of the Hamiltonian – one that flips the spin of the interacting fermions and the other that does not.

The 2-particle states consist of a $`q\overline{q}`$ pair sitting on the same $``$ lattice site. The 3-particle states consist of a $`q`$ and an $`\overline{q}`$ occupying neighboring sites separated by a link field (the gauge degrees of freedom in this approach). The 2-sector part of the Hamiltonian contains two terms that describe the kinetic energies of the $`q`$ and the $`\overline{q}`$, respectively. It also contains a longitudinally confining Coulomb interaction term acting between the pair. The 3-sector part of the Hamiltonian contains three terms that describe the kinetic energies of the $`q`$, the $`\overline{q}`$ and the link field, respectively. There are two more terms – one describes a longitudinally confining Coulomb interaction between the $`q`$ and the link field at the site where the $`q`$ resides and another that describes a longitudinally confining Coulomb interaction between the $`\overline{q}`$ and the link field at the neighboring site where the $`\overline{q}`$ resides. Note that there are no Coulomb terms connecting the $`q`$ and the $`\overline{q}`$ on neighboring sites. This is because these interaction terms are local on the $``$ lattice sites. Further, such interactions would yield non-planar diagrams which would violate our large $`N_c`$ approximation.

Transverse propagation of the mesons take place through the two interaction terms that connect the two Fock states. The meson propagates by emitting (2- to 3-sector) or absorbing (3- to 2-sector) a link field. When a meson is in the 2-sector (i.e., $`q`$ and $`\overline{q}`$ are sitting on the same $``$ lattice site), one of it’s constituents can hop to a neighboring site emitting a link field either through the spin flip term in which case it’s spin will flip or through the spin non-flip term in which case it’s spin will not flip. The meson in this configuration (a $`q`$ and an $`\overline{q}`$ on neighboring sites separated by one link field) represents it’s 3-particle state. When the meson is in this state, one of the fermion degrees of freedom hops and joins the other fermion on a neighboring site by absorbing the link field either through the spin flip in which case it’s spin will flip or through the spin non-flip interaction term in which case it’s spin will not flip. The nature of the propagation is reminiscent of that of an inchworm on a scale of a femtometer which motivated the name femtoworm.

The strength of spin flip and spin non-flip interaction terms are labeled $`m_v`$ and $`m_r`$, respectively. Though, in general, a link field emitted through a spin flip (spin non-flip) interaction term can be absorbed through a spin non-flip (spin flip) interaction term, such an interference term is disallowed in the $`𝐤_{}=\mathrm{𝟎}`$ case. This and other spin dynamics are discussed in details in .

The parameter space that defines the above $`q\overline{q}`$ bound state Hamiltonian contain the following parameters:
$``$ kinetic masses for the light quark in the 2 and 3 particle Fock sector.
$``$ kinetic masses for the strange quark in the 2 and 3 particle Fock sector.
$``$ kinetic mass for the link field.
$``$ spin flip and non-flip hopping terms for the light quark.
$``$ spin flip and non-flip hopping terms for the strange quark and
$``$ longitudinal gauge coupling.

As explained in , this parameter space is constrained by:
$``$ imposing Lorentz symmetry. This is done by demanding that the eigenstates of the Hamiltonian have the same dispersion relation.
$``$ imposing chiral symmetry in the spontaneously broken form. This is implemented by tuning the spin flip and spin non-slip hopping terms to yield the observed $`\pi \rho `$ mass splitting and
$``$ fixing the physical string tension. We used the string tension as obtained from glueball calculations .

## 3 Results and Conclusions

We have studied the LF momentum distribution amplitudes of several light mesons ($`\pi `$, $`\rho `$), mesons with strange quarks ($`KK^{}`$ system, $`\varphi `$) as well as heavy mesons ($`B`$-meson). For each of the above, we have also calculated their respective parton distribution functions.

The LF momentum distribution of the pion \[Fig.1(a)\] is quite close to the asymptotic shape. The first five even coefficients of a Gegenbauer projection of the momentum distribution, $`\psi (x)=4x(1x)_{i=0}^{\mathrm{}}a_iG_i^{\frac{3}{2}}(12x)`$, are $`\{1.5,\mathrm{\hspace{0.33em}0.102},0.006,0.005,0.001\}`$. The tail of the parton distribution function \[Fig.1(b)\] is somewhat higher than what is observed. Because there is only one link, most of the longitudinal momenta close to $`x=1`$ are carried by the fermions.

The dashed line in Fig.1(c) shows the $`\sqrt{z}`$ behavior of the $`B`$-meson LF momentum distribution in the heavy quark limit. For phenomenological applications  it is useful to have numerical estimates for the moments  (the normalization is $`_0^{\mathrm{}}𝑑z\varphi _{\mathrm{}}(z)𝑑z=1`$) for which we find $`𝑑zz\varphi _{\mathrm{}}(z)=1.51GeV`$ and $`𝑑zz^1\varphi _{\mathrm{}}(z)=1.22(GeV)^1.`$ For the decay constant , which plays also an important role in mixing phenomenology, we find $`f_B240MeV\pm 20MeV`$. Fig. 2 shows the Isgur-Wise form factor of the $`B`$-meson in the heavy quark limit which allows one to extract the $`V_{cb}`$ element of the CKM matrix from the zero recoil point of a semi-leptonic decay channels like $`B\overline{D}^{}l\nu `$. We obtain values of $`V_{cb}`$ in the range $`0.04150.0420`$ which is very consistent with the experimental data .

Since we used basis function techniques and extrapolated the results in Hilbert space dimensions, uncertainties from Hilbert space truncation are under control. Systematic errors arising from the Fock space truncation cannot be estimated unless higher Fock components are included. Thus, a main improvement over this work would be to include such Fock components. For example, we expect that including such components will significantly reduce the 2-particle Fock space amplitude, thereby reducing the decay constants which are calculated  from the normalizations of the LF distribution amplitudes. Other observables are expected to be effected in a similar manner.

For the heavy quark observables, $`1/M_b`$ and perturbative corrections are expected to contribute at zero recoil ($`v.v^{}=0`$). We were not able to estimate the former because we worked in the $`M_b\mathrm{}`$ limit. Also, it is not known how to reach small lattice spacings on the $``$ lattice and as such an estimate of the latter corrections are beyond the scope of this technique. Despite this caveat, we feel that with further improvements, this technique holds promise.


text:
# Spatial dispersion in CaF2 caused by the vicinity of an excitonic bound state.

## I INTRODUCTION

There is no doubt that CaF2 is the key component for establishing 157nm UV lithography. Its perfect cubic symmetry together with a great homogeneity of large single crystals which are very stable under normal conditions makes it a perfect material for homogeneous optical lens design.

However, there is a severe shortcoming of the material on which the attention was focused within the last year. This is the so called “intrinsic birefringence”.burnett2001  For the first glance it is astonishing. Looking at standard textbooks (e.g nye79 ) we note that a cubic crystal does not show any birefringence. On the other hand experiments observe a spatial dependence of the dielectric properties. Light propagating along different directions through the crystal shows a deviation between different polarization directions which is of the order of 10-6 nm/cm for UV light with 157nm. Or with other words the dielectric function depends on the direction in which the light propagates through the crystal and shows a difference which is of the order of 10-6.

How can this obvious misfit between standard textbooks and recent observations in CaF2 be resolved? The answer leads back to a physical phenomenon which is known as spatial dispersion.landaulifshitz8

The statement “a cubic system does never show birefringence” can be proven exactly, since in a cubic system the ellipsoid corresponding to the dielectric function is always a sphere. The question remains: What is a “system”? The crystal of CaF2 has exactly cubic symmetry but the “system” of light interacting with this perfectly cubic crystal shows a deviation from cubic symmetry. This is due to the fact that the light carries a $`𝐪`$ vector pointing in the vacuum in the direction of propagation. This is a consideration which has been already made more than one century ago lorentz1875 . However the speed of light is much larger than any speed in a solid state system and therefore one can usually safely neglect the wave-vector dependence of the propagating light. This statement is identical to the wavelength of light being in the optical window much larger than any length scale in a solid.

There are however two important exceptions where the assumption made above ceases to be valid. The first (I) is when one is by far leaving the optical window. This is the case at X-ray scattering experiments. Here the wavelength of light is of the order of the lattice constant and this causes a wave-vector (spatial) dependence of the resulting scattering signal. The second case (II) which is valid for CaF2 and many other solids is the case of light in the vicinity of a dispersive absorption occurring also at small wavevectors (large length scales). Since the dielectric function is an analytic function it has to follow causality which is expressed in the Kramers-Kronig relation between its real and imaginary part. The vicinity of a band edge can cause spatial dispersion in a semiconductor system (see e.g. cardona71 ).

In a semiconductor, the absorption is mainly dominated by the band edge and an exciton, the two particle bound state of the electron–hole pair only plays a minor role. The situation changes if we look at ionic crystals like CaF2. Here screening of the bare coulomb interaction is less effective which expresses itself in deep excitonic bound states. This is expressed in the exciton intensity being inverse proportional to the third power of the static dielectric function $`Iϵ_0^3`$.phillips70  The connection between a dispersive excitonic bound state and spatial dispersion has been already discussed for KI. meseguer84

When discussing possible alternative material designs to CaF2 it is important to understand the microscopic mechanism beyond the spatial dispersion in great detail and it is in particular important to know and predict the exciton dispersion.

The paper is organized as follows: In section II we define the system in the case when the approximation of locality breaks down. In section III we review the connection between the electronic part of the dielectric function and the bare two-particle propagator. We look at the free two particle propagator in the absence of any screened coulomb attraction and its connection to the band structure. The connection between the band structure, which can e.g. be based on DFT (density functional theory) calculations, and the dynamic polarization is investigated in section IV. In the following (section V) we discuss analyticity based on the resulting pole structure of the dielectric function. This describes to a large part the situation in a semiconductor where no deep bound state is present. Section (VI) serves to introduce the two body problem of the electron hole pair which leads to the excitonic bound state when a reasonable assumption for the weakly screened coulomb attraction has been made. Based on the model assumptions made, we discuss the spatial dispersion of CaF2 close to the bound state level in section VII and compare our predictions with the experimental work of.burnett2001

## II System of CaF2 interacting with light

Light has a $`𝐪`$ vector pointing along its propagation direction. This wave-vector becomes important when the light is interacting with a long wavelength ($`𝐪0`$) excitation of the system. In an isotropic medium this wave vector dependence will only depend on the absolute value of $`𝐪`$. In liquids and glasses this leads e.g. to the observation of hydrodynamic modes (phonon modes) in light scattering experiments.latz01  In single crystals the situation is different as not only the dependence on the absolute value of the wavevector of the light is important but also the direction of the wave vector with respect to the crystal orientations. The crystal of CaF2 possesses ideal cubic symmetry. That means in the static limit that also all response functions have to show this symmetry. Dynamically there are situations especially close to a strong excitation at long wavelengths, when the breaking of cubic symmetry due to the incident light gets important. Therefore not only the dynamic of a dielectric response has to be considered but also the wave-vector dependence. The electric field $`𝐄`$ is connected with the electric induction $`𝐃`$ in a non-local way.

$$𝐃(𝐪,\omega )=ϵ(𝐪,\omega )𝐄(𝐪,\omega )$$
(1)
where $`ϵ(𝐪,\omega )`$ is the dynamic dielectric function. For the physics of the spatial dispersion in CaF2 we will focus on the part of the dynamic dielectric function which is caused by electronic excitations. In general the dielectric function is an analytic function defined on the complex $`\omega `$-plane. Along the real axis it has an real and imaginary part.

$$ϵ(𝐪,\omega )=ϵ^{}(𝐪,\omega )+iϵ^{\prime \prime }(𝐪,\omega )$$
(2)
Only in the static limit $`\omega 0`$ one can assign a physical value to the real and imaginary part separately. The real part is connected with polarizabilities and the imaginary part with mobile charges ($`e^{\prime \prime }(\omega )=lim_{\omega 0}4\pi \sigma (\omega )/\omega `$) which are absent in an insulator. For a general crystal the dielectric function is a tensor. Since we want to keep the equations as transparent as possible we do not write tensor indices. We assume a coordinate system placed along the main axes of the cubic crystal. For a cubic crystal like CaF2 we can concentrate on specific directions of the $`𝐪`$ vector. Symmetry considerations allow to decompose all tensorial quantities into its symmetry components. This will provide rules to fill all the components of the tensorial dielectric function.

## III The dynamic dielectric function

The optical frequency regime is above the infrared optical phonon modes and below excitations over the band gap. In this regime the dynamical dielectric function can be decomposed in two parts.

$$ϵ(𝐪,\omega )=ϵ_{\mathrm{}}+ϵ_{e^{}h^+}(𝐪,\omega )$$
(3)
Where $`ϵ_{\mathrm{}}`$ is the IR-dielectric function which in a cubic crystal like CaF2 will be not only real but completely isotropic. The subscript $`\mathrm{}`$ refers to frequencies above the IR optical phonon modes but is still far below exciatations over the bandgap. Wave-vector dependence only enters from the second part which stems from electron hole excitations. The electron-hole dielectric function is proportional to the proper electron-hole polarization propagator fetwal71

$$ϵ_{e^{}h^+}(𝐪,\omega )\mathrm{\Pi }^{}(𝐪,\omega )$$
(4)
Where $`\mathrm{\Pi }^{}(𝐪,\omega )`$ is the propagator of an electron hole pair with total momentum $`𝐪`$. In the following we do not want to focus on details but rather want to outline the pole structure of this two particle propagator since it gives further insight into the problem.

## IV Connection to the band structure

If we look for the first glance at the bare electron hole propagator $`\mathrm{\Pi }^0(𝐪,\omega )`$ we can make the connection to the band structure.

$$\mathrm{\Pi }^0(𝐪,\omega )=\frac{1}{N}\underset{𝐤}{}\frac{f(𝐤,𝐪,\omega )}{\omega (E^e^{}(𝐤+𝐪)E^{h^+}(𝐤))}$$
(5)
The detailed form of the numerator can be taken from textbooks fetwal71  while the pole structure is connected with the electronic band structure. In general Eq. (5) should also contain a sum over all bands. Since we are approaching the band gap from the low frequency site and we know that CaF2 has a direct band gap we are for a first glance only interested in the highest valence band $`E^{h^+}(𝐤)`$ and the lowest conduction band $`E^e^{}(𝐤)`$ around the $`\mathrm{\Gamma }`$ point at $`𝐤=(0,0,0)`$. Even in a cubic crystal the band structure will show a different $`𝐤`$ dependence when looking along different directions in the crystal. In particular the $`(1,0,0)`$ direction will differ from the $`(1,1,0)`$ and $`(1,1,1)`$ direction.

## V Analytic properties of the dielectric function

From the analyticity of $`\mathrm{\Pi }^0(𝐪,\omega )`$ we can derive the qualitative behaviour of its real and imaginary part. This is plotted schematically in Fig. 1. In Fig. 1a the structure is shown in the $`\omega q`$ plane. There is the large band gap. An excitation over the bandgap is the creation of an electron-hole pair. Above the band gap there is a continuum of one particle excitations. The light dispersion is shown as well. Since the speed of light is much larger than any electronic speed one can usually in the optical frequency range safely neglect the wave-vector dependence of the light and use $`𝐪=0`$. This breaks down when excitations are considered which come close to the band gap. In this case it can become important that the band curvature is different along different directions in the crystal. Schematically we have plotted here a ficticious $`(1,1,0)`$ and $`(1,0,0)`$ band curvature.

In Fig. 1b the resulting real and imaginary part of $`\mathrm{\Pi }^0(𝐪,\omega )`$ is plotted resulting from the Kramers-Kronig relations. In a three dimensional system the imaginary part of the one particle continuum shows a $`\sqrt{\omega \omega _0}`$ behaviour. Therefore the real part does not diverge at the band edge (opposite to a 1D or 2D system) but increases up to a finite value. In order to keep the plots in Fig. 1 as simple as possible we have not plotted the real part of the polarizability above the gap.

If for different wave-vector directions in the crystal the band gap is shifted one can observe a spatial dispersion. This is the case in semiconductors where no deep excitonic bound state exists.cardona71  Here the spatial dispersion is caused by the wave-vector dependence of the one particle continuum. A weak attractive electron-hole interaction will not change the picture as long as no strong excitonic bound states are created below the one particle continuum. Therefore in these cases the dielectric function can already be described using $`ϵ\mathrm{\Pi }_0`$. The situation is completely different in a strongly ionic crystal like CaF2 here the excitonic bound state becomes important.

## VI The two particle bound state of the exciton

In a strongly ionic crystal screening is less effective and the coulomb attraction between the positively charged hole and the negatively charged electron leads to the formation of a strong two particle bound state. The solution of this two body problem involves in general a solution of the Bethe-Salpeter equation (see e.g. rohlfing98 ). However the basic principal of this bound state formation can already be understood with a much simpler interaction. The full solution for a complicated electron hole interaction for CaF2 has been taken into account in benedict99  here we just want to highlight the basic principles under which such a solution leads to an excitonic bound state. When we assume the attractive interaction to be only a function of the total momentum of the electron-hole pair $`U(𝐪)`$, the equations can be solved easily and one gets a RPA (random phase approximation) denominator for the proper polarization which stems from a summation of a geometrical series.

$$\mathrm{\Pi }^{}(𝐪,\omega )=\frac{\mathrm{\Pi }^0(𝐪,\omega )}{1U(𝐪)\mathrm{\Pi }^0(𝐪,\omega )}$$
(6)
This is also known as the ladder approximation to the Bethe-Salpeter equation. New physics will occur when the denominator of Eq. (6) shows additional poles beside the one particle continuum. This is the case when the denominator of Eq. (6) gets zero.

$$1U(𝐪)\mathrm{\Pi }^0(𝐪,\omega )=0$$
(7)
If this happens outside the one particle continuum a true bound state, the exciton is formed. Finding the bound state level can be schematically seen in Fig. 1b. When drawing a line at $`f(\omega )=1/U(𝐪)`$ it can for large $`U`$, for large attractive interactions cross the real part of $`\mathrm{\Pi }^0`$. This crossing marks the bound state level where a new pole occurs. Such a situation with a bound state is schematically sketched in Fig. 2. Below the one particle continuum the exciton level occurs as a true pole. The dynamic polarization and therefore also the dielectric function has to follow Kramers-Kronig relations connecting its real and imaginary part. In the vicinity of a true pole the real part of the dielectric function has to follow the asymptotic behaviour as

$$ϵ_{e^{}h^+}^{}(𝐪,\omega )\mathrm{}\left(\mathrm{\Pi }^{}(𝐪,\omega )\right)\frac{1}{\omega \omega _{ex}(𝐪)}=\frac{2\pi c}{\frac{1}{\lambda }\frac{1}{\lambda _{ex}(𝐪)}}=\frac{2\pi c\lambda \lambda _{ex}(𝐪)}{\lambda _{ex}(𝐪)\lambda }\frac{1}{\lambda \lambda _{ex}(𝐪)}$$
(8)
The wave-vector dependence of the exciton level $`\omega _{ex}(𝐪)`$ (or $`\lambda _{ex}(𝐪)`$) will in general differ along different directions in the crystal. This is schematically plotted in Fig. 2a. Exactly this difference causes the spatial dispersion seen in ionic crystals like CaF2 and BaF2.

## VII Spatial dispersion in CaF2 and predictions for experimental observations

In view of the measurements performed by Burnett et al. burnett2001  we need the refractive index along different directions. In an isotropic medium this is where the dielectric function is diagonal with three identical diagonal elements this would be simply given by:

$$n(\omega )=\sqrt{ϵ(\omega )}$$
(9)
with $`ϵ`$ being the (real part of the) diagonal part of the diagonalized dielectric function.

The concept of a refractive index becomes by far more sophisticated as soon as the medium shows anisotropies. This is also the case for CaF2. Here one has to go into a new coordinate system whose z-direction is determined by the $`𝐪`$-direction of the incident light. Two possible directions perpendicular to the z-direction complete this new coordinate system and are denoted with Greek letters $`\alpha ,\beta `$. Within the $`\alpha `$-$`\beta `$-plane lies the polarization vector $`𝐞`$ pointing in the direction of the electric induction $`𝐃`$ of the light. Within this new coordinate system the calculation of the refractive index is given by the zero of a two dimensional determinant.landaulifshitz8

$$\text{Det}\left|n^2\delta _{\alpha ,\beta }ϵ_{\alpha ,\beta }^1(𝐪,\omega )\right|=0$$
(10)
In general the eigenvalue problem of Eq. (10) has two solutions which correspond to the main polarization axes.

There is a further difference to an anisotropic medium. This is the fact that the wave-vector $`𝐪`$ is not in general parallel to the Poynting-vector $`𝐒`$ which points along the direction of the group velocity which is the analogon of the path in classical geometrical optic. This is expected to be small in a cubic crystal.

So the refractive index is dependent on $`𝐪`$ and can for a given experimental setup be dependent on the polarization direction of the incident light.

$$n=n(𝐪,\omega )=n(𝐪,2\pi c/\lambda )$$
(11)
In a cubic crystal like CaF2 there are several exceptions for the wave-vector $`𝐪`$ and the Poynting vector pointing along the same direction where no dependence on the polarization direction occurs as one can derive from symmetry considerations. One is for example light propagating along the main axis of the crystal (e.g. $`𝐪=(1,0,0)`$ now measured in the coordinate system of the crystal axes).

The dielectric function diverges close to the excitonic bound state depending on the wave-vector $`𝐪`$ according to Eq. (8) with $`\lambda _{ex}(𝐪)`$ is the wavelength of the exciton for the given $`𝐪`$ vector. This wavelength must be for CaF2 close to the 112nm where the exciton is found experimentally tomiki69 .

Note that Eq. (8) is quite different to what would result from just a formal Taylor expansion of the dielectric function. In this case the dielectric function should show a $`1/\lambda ^2`$ behaviour.

We want to look at the quantity measured in burnett2001 . There the difference of the refractive index $`\mathrm{\Delta }n`$ with a wave-vector pointing in the $`(1,1,0)`$ direction and the polarization vector $`𝐞`$ pointing in two different directions (in $`(1,1,0)`$ and $`(0,0,1)`$) is measured. Using Eq. (9) gives

$`2n\mathrm{\Delta }n`$ $``$ $`\mathrm{\Delta }ϵ`$ (12)
$`\mathrm{\Delta }n`$ $``$ $`{\displaystyle \frac{\mathrm{\Delta }ϵ}{2n}}{\displaystyle \frac{1}{\lambda \lambda _{ex(𝐪)}}}`$ (13)
In Fig. 3 we have plotted the wavelength dependence of $`\mathrm{\Delta }n`$ as it would result from the experimentally observed bound state level lying at 112nm.tomiki69

It can be seen that describing the experimental observations with Eq. (13) (dotted line) gives a better agreement with the experimental observation (crosses) than just a $`1/\lambda ^2`$ behaviour (dashed line) as it would be expected from a Taylor expansion of the dielectric function. We therefore claim in this work that the anisotropy which is seen in CaF2 at low wavelengths is caused by the exciton dispersion of an excitonic level around 112nm (11.2eV). Especially measurements in the wavelength range between 112nm and 157nm we predict to show a clear $`1/(\lambda \lambda _{ex}(𝐪))`$ behaviour.

## VIII conclusion

We argue that the optical anisotropy observed in CaF2 is caused by exciton dispersion. The delocalized exciton lies at approximately 112nm. Using this excitonic model the measured data on the optical anisotropy can be well explained and it is predicted that the dielectric function should diverge for lower wavelengths as $`1/(\lambda \lambda _{ex}(𝐪))`$. Here $`\lambda _{ex}(𝐪)`$ is the wavelength of the excitonic bound state.

If the connection to the exciton is indeed the dominant mechanism beyond the optical anisotropy all possible compensation approaches have to be reviewed. E.g. forming solid solution crystals, if possible at all, would only lead to an improvement, if the exciton dispersion will be smaller in such a crystal. Further the subtle difference between the spatial dispersion and birefringence has to be taken into account in a very careful way when performing lens design.

###### Acknowledgements.

M.L. thanks D.Strauch and K.Schmalzl for helpful discussions and for band structure data on CaF2. The work was supported by the BMBF project “Laserbasierte Ultrapräzisionstechnik – 157nm Lithographie, Teilvorhaben: Optische Materialien und Komponenten für die 157nm Lithographie, AP 4210”.


text:
# Definition 1.

A remark on the symmetries of the Riemann curvature tensor

The tetralinear form $`R(a,b,c,d)`$ defined as $`R(a,b)c,d`$, where $`R`$ is a Riemann curvature tensor, satisfies two basic symmetry properties:

1. $`R(a,b,c,d)=R(b,a,c,d)=R(a,b,d,c)`$,
2. $`R(a,b,c,d)+R(b,c,a,d)+R(c,a,b,d)=0`$.

Out of these one can derive other identities, e.g. $`R(a,b,c,d)=R(c,d,a,b)`$. The purpose of this note is to show that these properties can be (efficiently) expressed as $`SL(2)`$-invariance.

###### Definition 1.

If $`V`$ is a vector space then let $`(V)(V^{})^4`$ be the subspace of the tetralinear forms satisfying the properties 1. and 2.

###### Lemma 1.

Let $`(W,\rho )`$ be a finite-dimensional representation of $`SL(2)`$. If $`wW`$ is invariant with respect to the subgroup of upper-triangular matrices $`SB(2)SL(2)`$, then it is $`SL(2)`$-invariant.

###### Proof.

The $`SB(2)`$ invariance means that $`w`$ is a highest weight vector of weight 0. ∎

###### Proposition 1.

Let $`U`$ be an auxiliary 2-dimensional vector space with a chosen area form $`\omega `$. Then, if $`V`$ is a finite-dimensional vector space, there is a natural isomorphism

$$(V)\left(^4(UV^{})\right)^{SL(U)}.$$
(1)
###### Proof.

Let $`u_1`$, $`u_2`$ be a basis of $`U`$ such that $`\omega (u_1,u_2)=1`$ (the isomorphism will be constructed using this choice, but due to the $`SL(U)`$-invariance, it will not depend on it). It gives us an isomorphism $`UV^{}(u_1V^{})(u_2V^{})`$, i.e.

$$^4(UV^{})\underset{i=0}{\overset{4}{}}^i(u_1V^{})^{4i}(u_2V^{}).$$
(2)
If we restrict the left hand side to the elements invariant w.r.t. the diagonal matrices in $`SL(2)`$ (acting by $`u_1\mu u_1`$, $`u_2\mu ^1u_2`$), we get $`^2(u_1V^{})^2(u_2V^{})`$; this space can be identified with the tetralinear forms on $`V`$ satisfying 1. Restricting further to the elements invariant w.r.t. $`u_1u_1`$, $`u_2u_2+\lambda u_1`$ we get the forms satisfying both 1. and 2. We have passed to the subspace of $`SB(2)`$-invariants; by Lemma 1 this is the same as $`SL(2)`$-invariants. ∎

Final remarks. As a simple example, the identity $`R(a,b,c,d)=R(c,d,a,b)`$ follows from invariance w.r.t. the matrix $`u_1u_2`$, $`u_2u_1`$. More importantly, there is a supergeometry behind the ad hoc construction of this note, making it more natural; it will (hopefully) appear elsewhere.

Pavol Ševera

Dept. of Theor. Physics
Comenius University
Bratislava, Slovakia