text:
# The Size and Shape of Voids in Three-Dimensional Galaxy Surveys

## 1 Introduction

When the positions of galaxies are mapped in redshift space, the nature of the large-scale structure seen in the maps is “frothy” or “bubbly”. Voids – regions in which there are few or no galaxies – fill much of the map. The galaxies exist mainly in thin filaments or sheets that lie between the voids. Comparison of the observed pattern of voids with that predicted by a given model of structure formation is potentially a very powerful test for acceptance or rejection of that model. Unfortunately, the differences between models are often subtle. Many models for structure formation predict a bubbly galaxy distribution. In the gravitational instability scenario, regions that are originally underdense expand faster than the Hubble flow (Fillmore & Goldreich, 1984; Bertschinger, 1985). As underdense regions expand, they become more nearly spherical (Fujimoto, 1983; Icke, 1984; Blaes, Goldreich, & Villumsen, 1990). Whether voids in the actual universe are viewed as isolated spherical structures or as a space-filling foam depends on the density level at which voids are defined (Gott, Melott, & Dickinson, 1986). Van de Weygaert & van Kampen (1993) did extensive studies of voids in constrained Gaussian fields, finding that voids become nearly spherical in their very underdense inner regions, while their boundary regions remain more irregular.

If structure grows via gravitational instability, then the size and shape of voids depends on the initial power spectrum $`P(k)`$ for the density fluctuations and on the density parameter $`\mathrm{\Omega }`$. In an open universe, void evolution stops when $`\mathrm{\Omega }0.5`$; thereafter, the void network simply expands along with the Hubble flow (Regös & Geller, 1991). At a given epoch, the mean void radius is proportional to the nonlinearity scale $`k_{\mathrm{nl}}`$ (Kauffmann & Melott, 1992). The full distribution of void sizes, however, depends on the shape of the spectrum $`P(k)`$ (Melott, 1987; Ryden & Melott, 1996).

If we know with absolute accuracy the position of galaxies in real space, we could use the spectrum of void sizes to place constraints on the initial power spectrum $`P(k)`$. However, measuring the distances to galaxies is difficult; it is much easier to measure their redshifts. Consequently, practical studies of the properties of voids must take place in redshift space rather than real space. If all galaxies smoothly followed the Hubble flow, with no peculiar velocities, and if the Hubble constant $`H_0`$ were truly constant with time, then the mapping between real space and redshift space would be linear. Generally, though, the Hubble constant changes with time, and structure that is isotropic in real space becomes distorted along the line of sight in redshift space, with the amount of distortion increasing with redshift $`z`$ (Alcock & Paczyński, 1979). Potentially, distortions in the shapes of voids in redshift space can provide an estimate of the deceleration parameter $`q_0`$ (Ryden, 1995). However, the distortions that result from cosmological effects become large only when $`z1`$. In the nearby universe, where $`z<1`$, the dominant contribution to distortions from redshift space comes from the peculiar velocity of galaxies. Ryden & Melott (1996) demonstrated, for instance, that in two-dimensional simulations with $`P(k)k^n`$, the characteristic void size is larger in redshift space than in real space. Redshift space distortions can also, in principle, be harnessed to measure $`\mathrm{\Omega }_m`$, the density in nonrelativistic matter (Melott et al., 1998). One goal of this paper is to extend the analysis of Ryden & Melott (1996) to three-dimensional simulations with more realistic CDM power spectra.

Analyzing the properties of voids first requires a definition of what a void is. Some of the statistics used to describe voids define a void as a region totally devoid of galaxies. The void probability function (VPF) is one such statistic; the VPF $`P_0(V)`$ is the probability that a randomly positioned sphere of volume $`V`$ contains no galaxies (White, 1979). However, the random dilution of a point process (selecting only a fraction of all the galaxies in the universe, for instance), although it leaves the correlation function unchanged, strongly affects VPF (Sheth, 1996). Thus, from a practical viewpoint, it makes more sense to define a void which is underdense with respect to the average number density of galaxies. Little & Weinberg (1994), for instance, use the underdense probability function (UPF). The UPF $`P_f(V)`$ can be defined as the probability that a randomly positioned sphere of volume $`V`$ contains a number density of galaxies less than or equal to $`f`$ times the average number density of galaxies in the entire survey. We will follow Little & Weinberg (1994) in setting a void threshold of $`f=0.2`$, matching the density contrast of the largest voids in the CfA redshift survey (Vogeley, Geller, & Huchra, 1991; Vogeley et al., 1994). One advantage of the UPF over the VPF is that it is relatively insensitive to the sparseness with which the galaxies are sampled. In this paper, we will explicitly compare the UPF and VPF for different galaxy sampling densities.

In addition to statistical measurements such as the VPF and UPF, algorithms exist that identify individual voids within a sample (Kauffmann & Fairall, 1991; Kauffmann & Melott, 1992; Ryden, 1995; Ryden & Melott, 1996; El-Ad & Piran, 1997; Aikio & Mähönen, 1998). In addition to providing a spectrum of void sizes, these void-detection algorithms also enable us to specify the location and shape of individual voids. In this paper, we first use an algorithm (an extension of the two-dimensional algorithm of Ryden ) which defines voids as non-overlapping underdense spheres. Although this algorithm is conceptually simple, it has the disadvantage of forcing voids to be spherical. Thus, we also employ a more sophisticated algorithm, based on that of Aikio & Mähönen (1998; AM), which does not constrain voids to be any particular shape. The original AM algorithm defines voids as empty regions; our modification of the original defines voids as underdense regions. For the voids found by the modified AM algorithm, we introduce an elongation statistic $`Q`$ which measures whether the voids are stretched or squashed along the line of sight from the origin. Just as with the UPF and VPF, we also apply the void-detection algorithms to different sampling densities of the same survey. Our studies of void properties in numerical simulations will permit a more effective interpretation and understanding of voids in future three-dimensional redshift surveys.

In section 2, we describe the N-body simulation analyzed in this paper. In section 3, we examine the statistical properties of voids, as given by the VPF and UPF. In section 4, we apply void-detection algorithms to the simulations. Finally, in section 5, we analyze the effects of peculiar velocities and galaxy sampling density on void properties, and discuss implications for future work.

## 2 The Numerical Simulation

The simulation used for testing in this paper was done using a particle-mesh (PM) N-body simulation. The PM method is quite fast, and with a mean particle density of one per simulation cell, represents the maximum resolution that can be achieved without introducing two-body scattering that decouples the result from its initial conditions on small scales (Kuhlman, Melott, & Shandarin, 1996; Splinter et al., 1998). The simulation used here had 2563 particles on a 2563 mesh. Initial conditions were generated by Fast Fourier Transform with random phase perturbations. Since behavior at high density peaks is not of interest to us in our study of voids, the simulations were begun with an RMS density of 0.25 at the resolution limit. The simulations were evolved until the RMS overdensity inside a randomly located sphere of radius $`8h^1\mathrm{Mpc}`$ was $`\sigma _8=1.05`$.

A matter-dominated Friedman-Robertson-Walker background density was assumed, with the cosmological constant $`\mathrm{\Lambda }`$ set equal to zero. Since nonlinear modes are filtered out, and the dynamical effect of nonzero $`\mathrm{\Omega }_\mathrm{\Lambda }`$ is well-understood in perturbation theory (Lahav et al., 1991), we did not use nonzero values of $`\mathrm{\Lambda }`$. The value of $`\mathrm{\Omega }_m`$, the dimensionless matter density, was taken to be $`\mathrm{\Omega }_m=1`$. The box size was take to be 1536 Mpc and the Hubble constant $`h=H_0/100\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1=0.75`$. Thus, the box size in redshift space corresponds to $`1.152\times 10^5\mathrm{km}\mathrm{s}^1=0.384c`$. Formally, in the simulation $`h`$ only sets an overall timescale; since both the expansion rate and the particle velocities scale with $`h`$, the redshift space appearance does not change with $`h`$, only its overall scale.

The initial power spectrum we assume is a cold dark matter (CDM) spectrum, in which the parameter $`\mathrm{\Gamma }`$ determines the shape of the spectrum. Smaller values of $`\mathrm{\Gamma }`$ are favored today, because they push the turnover in the slope of the power spectrum to large scales, in better agreement with observations. To test our void-finding algorithms, we use $`\mathrm{\Gamma }=0.25`$. Normally, $`\mathrm{\Gamma }`$ is taken to be $`\mathrm{\Omega }_0h`$, as the turnover scale is set by the horizon at the end of the radiation-dominated era. However, we break this assumed coupling and take $`\mathrm{\Gamma }`$ as a free parameter descriptive of the spectral shape. Thus, our model corresponds to what is sometimes called $`\tau `$CDM.

The number density of mass points in the simulation, $`n=4.63\times 10^3\mathrm{Mpc}^3`$, corresponds to the number density of galaxies with $`L>0.42L_{}`$, assuming that the luminosity function of galaxies is a Schechter function with slope $`\alpha =1.07`$ and normalization $`\mathrm{\Phi }_{}=6.75\times 10^4\mathrm{Mpc}^3`$ (Efstathiou, Ellis, & Peterson, 1988). One purpose of this paper is to investigate the dependence of void properties on the sampling density of a galaxy survey. Thus, from the initial numerical simulation, we have created three different samples, corresponding to volume-limited surveys with different luminosity cutoffs. To create our densest sample, we randomly selected a fraction $`X=0.82`$ of the initial mass points, in order to match the number density of galaxies with $`L>L_{}/2`$. For the next densest sample, we selected a fraction $`X=0.31`$ of the mass points, to match the number density of galaxies with $`L>L_{}`$. Finally, for our least dense sample, we selected a fraction $`X=0.067`$ of the mass points, to match the density of galaxies with $`L>2L_{}`$. Note that by randomly selecting galaxies in this way, we are assuming that galaxies are unbiased with respect to the mass distribution, and that luminosity segregation does not exist. We wish to study only the effects of different sampling densities, and not the more subtle effects of bias.

For each sampling density, we create two mock galaxy surveys, one without peculiar velocities (which we will call the “real space” survey) and one with peculiar velocities (the “redshift space” survey). The mock surveys are spheres with a radius of 480 Mpc, or $`z=0.12`$ in redshift units. Note that in a flat $`\mathrm{\Omega }_m=1`$ universe, a galaxy with $`L=L_{}/2`$ located at $`z=0.12`$ will have an apparent magnitude of $`m_B19`$, assuming $`M_B^{}=19.7+5\mathrm{log}h`$ (Efstathiou, Ellis, & Peterson, 1988). For a typical galaxy color of $`r^{}B1`$, this roughly corresponds to the flux limit $`r^{}=17.7`$ expected for the galaxy redshift sample of the Sloan Digital Sky Survey (Weinberg, 2000). A slice through the $`L>L_{}`$ “real space” survey is shown in Figure 1a; a slice through the $`L>L_{}`$ “redshift space” survey is shown in Figure 1b.

## 3 Statistics of Voids

The void probability function, or VPF (White, 1979) has been a widely used statistic for measuring the characteristic size of voids. The VPF $`P_0(V)`$ is defined as the probability that a randomly located sphere of volume $`V`$ contains no galaxies. The VPF has been applied to numerical simulations (Fry et al., 1989; Einasto et al., 1991; Weinberg & Cole, 1992; Little & Weinberg, 1994; Vogeley et al., 1994; Ghigna et al., 1994, 1996; Colombi, Bouchet, & Hernquist, 1996; Kauffmann, Nusser, & Steinmetz, 1997) and to redshift surveys (Fry et al., 1989; Einasto et al., 1991; Vogeley, Geller, & Huchra, 1991; Vogeley et al., 1994; Ghigna et al., 1996). The VPF for the simulation studied in this paper is given in the left-hand panels of Figure 2. The volumes plotted are in redshift units; to convert to physical units, multiply by $`(c/H_0)^3=6.4\times 10^{10}\mathrm{Mpc}^3`$. In each panel, the solid line shows the VPF in real space, and the dotted line shows the VPF in redshift space. The effect of peculiar velocity distortions is to increase the void probability at a given volume $`V`$. A comparison of the VPF for the different sampling densities, however, vividly illustrates the very strong dependence of the VPF on the mean interparticle spacing. To illustrate, let’s define a characteristic void volume $`V_{\mathrm{VPF}}`$ as the volume for which $`P_0(V_{\mathrm{VPF}})=0.01`$. For a Poisson distribution of points, $`P_0(V)=\mathrm{exp}(\nu V)`$, where $`\nu `$ is the mean number density, and thus $`V_{\mathrm{VPF}}=4.61/\nu `$. For our $`L>L_{}/2`$ sample, $`V_{\mathrm{VPF}}=4.0\times 10^8=9.6/\nu `$; for our $`L>L_{}`$ sample, $`V_{\mathrm{VPF}}=1.14\times 10^7=10.3/\nu `$; and for our sparse $`L>2L_{}`$ sample, $`V_{\mathrm{VPF}}=4.3\times 10^7=8.3/\nu `$. Thus, the voids in the simulation are larger than those in a Poisson distribution of equal $`\nu `$, but it is still approximately true that the characteristic void size (defining voids as totally empty volumes) is proportional to $`\nu ^1`$.

The dependence of characteristic void size on the sample density $`\nu `$ is reduced if we use the underdense probability function (UPF) to measure the statistics of voids. The UPF $`P_f(V)`$ is defined as the probability that a randomly located sphere of volume $`V`$ contains a number density of galaxies that is less than $`f\nu `$, where $`0<f<1`$, and $`\nu `$ is the mean number density of galaxies in the sample. (In a flux-limited survey, we can generalize the UPF so that $`P_f(V)`$ is the probability that a sphere of volume $`V`$ located at a distance $`r`$ from the origin contains a number density of galaxies that is less than $`f\nu (r)`$, where $`\nu (r)`$ is the mean number density of detected galaxies at $`r`$.) Following previous work (Vogeley, Geller, & Huchra, 1991; Weinberg & Cole, 1992; Ryden & Melott, 1996), we will set our density threshold at $`f=0.2`$. The UPF for the simulation studied in this paper is given in the right-hand column of Figure 2; in each panel, the solid line shows the UPF in real space and the dotted line shows the UPF in redshift space. Again, the effect of peculiar velocity distortions is to increase the probability of finding a void of volume $`V`$. For volumes $`V<1/(f\nu )`$, the UPF and the VPF are identical, since such a small volume can only fall below the density threshold if it contains no galaxies. In general, a sphere of volume $`V`$ will be underdense if it contains at most $`M`$ galaxies, where $`M=\mathrm{int}(f\nu V)`$. When $`M`$ is small, the UPF shows discreteness effects, visible as the sawtooth pattern in Figure 2; the UPF jumps upward whenever $`V`$ is an integral multiple of $`1/(f\nu )`$.

The characteristic void size $`V_{\mathrm{UPF}}`$ is defined as the volume for which $`P_{0.2}(V_{\mathrm{UPF}})=0.01`$. (The choice of $`P_0=0.01`$ as the defining probability is somewhat arbitrary; however, because of the rapid decline of $`P_0`$ at large volumes, the exact value of $`P_0`$ is not crucial.) For our dense $`L>L_{}/2`$ sample, $`V_{\mathrm{UPF}}=2.7\times 10^7=64/\nu `$; for the $`L>L_{}`$ sample, $`V_{\mathrm{UPF}}=4.2\times 10^7=38/\nu `$; for the sparse $`L>2L_{}`$ sample, $`V_{\mathrm{UPF}}=9.2\times 10^7=17.7/\nu `$. The inclusion of peculiar velocity distortions increases $`V_{\mathrm{UPF}}`$ by a factor which ranges from $`2.5`$ for the densest sample to $`1.4`$ for the sparsest sample. Although the dependence of the underdense void size $`V_{\mathrm{UPF}}`$ on $`\nu `$ is less strong than that of $`V_{\mathrm{VPF}}`$, it is not true that $`V_{\mathrm{UPF}}`$ is independent of $`\nu `$ for plausible sampling densities. Thus, in a flux-limited survey, where the measured $`\nu `$ decreases with distance from the origin, the characteristic void size $`V_{\mathrm{UPF}}`$ will increase with distance.

## 4 Void-detection Algorithms

The UPF gives a statistical measure of the number and size of voids in a given galaxy distribution. Frequently, however, it is useful to identify individual voids instead of simply giving a statistical description. Many different algorithms have been used to detect and identify individual voids (Kauffmann & Fairall, 1991; Kauffmann & Melott, 1992; Ryden, 1995; Ryden & Melott, 1996; El-Ad & Piran, 1997; Aikio & Mähönen, 1998). In this paper, we will investigate two void-detection algorithms that are distinguished by their ease of use and clarity of conception. The first algorithm, based on that of Ryden (1995), identifies voids as nonoverlapping spheres; we will call this the “sphere algorithm”. The second algorithm, based on that of Aikio and Mähon̈en (1998; AM) permits voids to be nonspherical; we will call this the “AM algorithm”.

Both void-detection algorithms start by defining a continuous scalar field $`D_f(\stackrel{}{x})`$ within the galaxy survey. At any location $`\stackrel{}{x}`$, $`D_f`$ is defined as the radius of the largest sphere centered on $`\stackrel{}{x}`$ within which the average galaxy density is equal to $`f\nu `$. To implement the “sphere algorithm”, first locate the global maximum of $`D_f`$ within the survey; call the location of this maximum $`\stackrel{}{x}_1`$. This is the center of the largest spherical void in the sample, which has a radius $`D_f(\stackrel{}{x}_1)`$. To find the second largest void, find the point $`\stackrel{}{x}_2`$ for which $`D_f`$ is maximized, subject to the constraint that

$$D_f(\stackrel{}{x}_2)+D_f(\stackrel{}{x}_1)|\stackrel{}{x}_2\stackrel{}{x}_1|.$$
(1)
The point $`\stackrel{}{x}_2`$ is then the center of the second-largest spherical void, which has radius $`D_f(\stackrel{}{x}_2)`$. In other words, the second-largest void is the largest underdense sphere that doesn’t overlap the largest void. Additional voids are found by an iterative process. The $`N^{\mathrm{th}}`$ largest void is located at the position $`\stackrel{}{x}_N`$ for which $`D_f`$ is maximized, subject to the constraint that

$$D_f(\stackrel{}{x}_N)+D_f(\stackrel{}{x}_i)|\stackrel{}{x}_N\stackrel{}{x}_i|,$$
(2)
for $`i=1,2,\mathrm{},N1`$.

In practice, we compute the values of $`D_f`$ on a cartesian grid superimposed on the galaxy distribution. For the simulation used in this paper, the grid spacing we used was $`\mathrm{\Delta }x=0.0012`$ in redshift units ($`\mathrm{\Delta }x=4.8\mathrm{Mpc}`$ in physical units). This is 0.8 times the resolution of the original numerical simulation. Using a grid very much finer than that of the original simulation is pointless, since there is no information on such small scales. The computed value of $`D_f`$ for each grid point was forbidden to be larger than the distance from the grid point to the sample boundary at z=0.12. We then located the spherical voids, using the algorithm outlined above, subject to the additional constraint that the void centers lie on grid points. Since the discreteness of the superimposed grid creates errors of order $`\mathrm{\Delta }x`$ in the location of void centers, we halt the void-detection algorithm when $`D_f=2\mathrm{\Delta }x`$.

To emphasize the difference between the properties of totally empty voids and underdense voids, we implemented the sphere algorithm twice, once with $`f=0`$, and once with $`f=0.2`$. To illustrate the voids found with the sphere algorithm, Figure 3 shows a slice through the underdense spheres found in the $`L>L_{}`$ sample, setting $`f=0.2`$. Figure 3a shows the spherical voids found in real space (without peculiar velocity distortions) and Figure 3b shows the spherical voids found in redshift space (with peculiar velocity distortions).

To show the spectrum of void sizes found with the sphere algorithm, Figure 4 plots $`F_f(V)`$, the fraction of the total volume of the sample found in spherical voids with volume $`V`$. The right column of Figure 4 displays $`F_{0.2}`$, the fraction of the total volume found in voids with underdensity $`f=0.2`$. For comparison, the left column of Figure 4 gives $`F_0`$, the fraction of the total volume found in totally empty voids. The solid line in each panel gives the distribution in real space, and the dotted line gives the distribution in redshift space. The left column of Figure 4 demonstrates, as in the case of the VPF, that the characteristic size of empty voids depends strongly on the density $`\nu `$ of the sample. Because of this undesired feature, we will only examine in detail the properties of underdense voids.

A characteristic void size can be found, following the practice of Kauffmann & Melott (1992), by computing the volume-weighted mean void size,

$$\overline{V}=\underset{1}{\overset{N_V}{}}V_i^2/\underset{1}{\overset{N_V}{}}V_i,$$
(3)
where $`V_i`$ is the volume of the $`i^{\mathrm{th}}`$ largest void, and the total number of voids is $`N_V`$. Alternatively, we can define a characteristic void size $`V_X`$ which is the void volume such that $`X`$ percent of the total volume in the sample is contained in voids of size $`V_X`$ or bigger. The characteristic void sizes $`\overline{V}`$ and $`V_{20}`$ are given in Table 1 for underdense voids in samples of different number density, with and without peculiar velocities.

The characteristic size of spherical voids (whether defined as $`\overline{V}`$ or $`V_{20}`$) increases in going from real space to redshift space. The more densely sampled the survey, the greater the increase in void size. For the $`L>2L_{}`$ sample, $`\overline{V}`$ increases by 43% in going from real space to redshift space, and $`V_{20}`$ increases by 40%. For the $`L>L_{}/2`$ sample, $`\overline{V}`$ increases by 88% in going from real space to redshift space, and $`V_{20}`$ by 155%.

Although the sphere algorithm gives a rough estimate of the spectrum of void sizes, the volumes of the voids found by this algorithm will generally be underestimates of the “true” void size. That is, if a sample contains an empty region surrounded by a continuous, well-defined, extremely overdense wall, the void found by the sphere algorithm will be the largest sphere that can be inscribed within the wall. The remaining empty space within the void wall will then be iteratively filled with smaller and smaller spheres. A more flexible algorithm – one which doesn’t impose the artificial constraint that voids are spherical – should give a more accurate measure of void size.

In addition to underestimating the size of voids, the sphere algorithm gives no hint of the true void shape. In two dimensions, Ryden (1995) and Ryden & Melott (1996) estimated the shapes of voids by fitting ellipses to the underdense region. Each void was then characterized by an axis ratio $`q`$ and a position angle $`\varphi `$ in addition to its area $`a`$ and position $`(x,y)`$. Extending this algorithm to three dimensions by fitting ellipsoids to the underdense regions becomes a computationally daunting task. Each ellipsoidal void must be characterized by two axis ratios and three Euler angles in addition to its volume $`V`$ and position $`(x,y,z)`$. The introduction of additional parameters makes the search through parameter space far more time-consuming. Thus, instead of approximating the shape of voids by fitting ellipsoids to them, we adopted the more flexible scheme of Aikio & Mahönën (1998; AM).

To implement the AM algorithm, we start with the field $`D_f`$ as defined on the cartesian grid that we have superimposed on the galaxy survey. If the grid spacing is $`\mathrm{\Delta }x`$, then each grid point can be thought of as being at the center of a cubical “elementary cell” of volume $`(\mathrm{\Delta }x)^3`$. We identify the local maxima of $`D_f`$ on the grid as being those points which have values of $`D_f`$ greater than that of their 26 closest neighbors (the 6 points at a distance $`\mathrm{\Delta }x`$, the 12 points at a distance $`\sqrt{2}\mathrm{\Delta }x`$, and the 8 points at a distance $`\sqrt{3}\mathrm{\Delta }x`$). We label the local maxima we find as $`M_1`$, $`M_2`$, $`\mathrm{}`$, $`M_N`$, where $`N`$ is the total number of maxima located. The AM algorithm assigns every elementary cell to a “subvoid” associated with some maximum $`M_i`$. To discover which subvoid a particular elementary cell belongs to, a “climbing algorithm” is used. For a elementary cell $`X`$, we compute the gradient in $`D_f`$ to each of the neighboring cells. We then “climb” to elementary cell $`X^{}`$ for which the gradient has the largest (positive) value. The climbing continues from cell to cell until a local maximum $`M_i`$ is reached. The cell $`X`$ (and every other cell along the climbing route) is then assigned to the subvoid of maximum $`M_i`$. In this way, every elementary cell is assigned to a subvoid, and each maximum $`M_i`$ has an associated subvoid which consists of at least one elementary cell. Once every elementary cell is assigned to a subvoid, the subvoids are joined together into larger voids. The subvoids associated with maxima $`M_i`$ and $`M_j`$ are members of the same void if the distance between $`M_i`$ and $`M_j`$ is less than both $`D_f(M_i)`$ and $`D_f(M_j)`$. Using this criterion, all the subvoids are grouped into voids; some voids contain a single subvoid, while others contain many subvoids linked together in a friends-of-friends percolation. The AM algorithm and its implementation is described in more detail in the original paper by Aikio & Mähönen (1998); our modification is to use the underdensity field $`D_f`$, where AM restricted themselves to the case $`f=0`$.

Figure 5a shows a slice through the underdense voids in the $`L>L_{}`$ sample, without the inclusion of peculiar velocity distortions. The AM algorithm was used with a density threshold $`f=0.2`$. In each panel of Figure 5, an elementary cell is colored white if it belongs to the same void as its 26 nearest neighboring cells; it is colored black if one or more of its neighbors belongs to a different void. Figure 5b shows a slice through the underdense voids in the $`L>L_{}`$ sample, this time with the effects of peculiar velocities included.

The volume $`V`$ of an individual void found by the AM algorithm can be found by simply adding together the volume $`(\mathrm{\Delta }x)^3`$ of the elementary cells which it contains. The full spectrum of void sizes in a particular galaxy sample is given by $`F_f(V)`$, the fraction of the total volume of the sample found in voids with volume $`V`$. The left column of Figure 6 shows $`F_0`$, the distribution of void volumes using the original algorithm of AM, in which $`f=0`$. The right column of Figure 6 displays $`F_{0.2}`$, using our usual underdense criterion, $`f=0.2`$. Some properties of the AM voids are the similar to those of the spherical voids found earlier. For instance, the left column of Figure 6 shows that the characteristic size of empty AM voids (just like those of spherical voids) is strongly dependent on the density $`\nu `$ of the sample. Also, peculiar velocities increase the characteristic size of AM voids as well as of spherical voids. Using the AM algorithm with $`f=0.2`$, the volume-weighted mean void size $`\overline{V}`$ and the characteristic size $`V_{20}`$ are given in Table 1. For the sparsely sampled survey ($`L>2L_{}`$), $`\overline{V}`$ increases by 58% in going from real space to redshift space, and $`V_{20}`$ increases by 66%. For the densely sampled survey ($`L>L_{}/2`$), $`\overline{V}`$ increases by 79% in going from real space to redshift space, and $`V_{20}`$ increases by 82%. In addition, AM voids as well as spherical voids are larger at lower sampling densities. Given a density threshold of $`f=0.2`$, the characteristic void size (either $`\overline{V}`$ or $`V_{20}`$) for the sparsely sampled survey ($`L>2L_{}`$) is 2.5 times that of the densely sampled survey ($`L>L_{}/2`$).

One important difference between the voids found by the AM algorithm and those found by the sphere algorithm is that the AM voids are larger. At all sampling densities, it is found that $`\overline{V}(\mathrm{AM})/\overline{V}(\mathrm{sphere})2.7`$. The ratio $`V_{20}(\mathrm{AM})/V_{20}(\mathrm{sphere})`$ ranges from $`6`$ for the $`L>2L_{}`$ survey through $`8`$ for the $`L>L_{}`$ survey to $`14`$ for the $`L>L_{}/2`$. Thus, the answer to the seemingly innocuous question, “How large is a typical void?” depends not only on the sampling density of the galaxy survey and on whether the survey is done in real space or in redshift space, but also on the void-finding algorithm used and on the definition adopted for the typical void size.

Another important difference between the voids found by the AM algorithm and those found by the sphere algorithm is that the AM voids are not compelled to be spherical. Hence the shapes of the AM voids can be used as a measure of the shape of voids in the galaxy distribution. It is of particular interest to discover whether voids are distorted along the line of sight from the observer at the origin. At relatively small redshifts ($`z1`$), the dominant source of distortion in redshift maps is the peculiar velocities of galaxies. In examining two-dimensional simulations with power spectra $`Pk^n`$, Ryden & Melott (1996) found only a mild tendency for voids to be distorted or compressed along the line of sight. With an $`n=2`$ power spectrum, the largest voids were slightly elongated along the line of sight; with an $`n=0`$ spectrum, voids were slightly compressed along the line of sight.

To find whether the AM voids are preferentially elongated or compressed along the line of sight, we begin by computing the moments of the voids. If a void contains $`N`$ elementary cells, with the center of the $`i^{\mathrm{th}}`$ cell at $`(x_i,y_i,z_i)`$, the moments of the void can be computed as

$$\overline{x^\alpha y^\beta z^\gamma }=\frac{1}{N}\underset{i=1}{\overset{N}{}}x_i^\alpha y_i^\beta z_i^\gamma .$$
(4)
The “center of mass” of the void, weighting all elementary cells equally, is at $`(\overline{x},\overline{y},\overline{z})`$. For a given void, we can create a new coordinate system, with its origin at $`(\overline{x},\overline{y},\overline{z})`$, with its $`x^{}`$ axis passing through the location of the observer, and with its $`y^{}`$ and $`z^{}`$ axis perpendicular to the $`x^{}`$ axis. If we know the coordinates $`(x_i,y_i,z_i)`$ of a mass element in the old coordinate system (centered on the observer), we can compute the coordinates $`(x_i^{},y_i^{},z_i^{})`$ in the new coordinate system (centered on the void center).

A dimensionless measure of a void’s elongation or compression along the line of sight is

$$Q\frac{3\overline{(x^{})^2}}{\overline{(x^{})^2}+\overline{(y^{})^2}+\overline{(z^{})^2}}.$$
(5)
The quantity $`Q`$ has some useful properties. Its denominator is independent of the orientation of the coordinate system; it’s simply the mean square distance of all the elementary cells from the void center. For a void of arbitrary shape, the mean value of $`Q`$, averaged over all viewing angles, is $`<Q>=1`$. Thus, for a population of voids oriented randomly with respect to the observer, we expect the average value of $`Q`$ to be 1. A value $`Q>1`$ indicates that a void is elongated along the line of sight; a value $`Q<1`$ indicates that a void is compressed along the line of sight. For any void, in any orientation, $`0Q3`$. As an example, consider a triaxial ellipsoid with principal axes of length $`abc`$. The maximum value of $`Q`$ for this ellipsoid occurs when the long axis is aligned with the line of sight from the observer to the ellipsoid’s center. In this case,

$$Q_{\mathrm{max}}=\frac{3a^2}{a^2+b^2+c^2}1.$$
(6)
The minimum value of $`Q`$ occurs when the short axis is aligned with the line of sight. In this case,

$$Q_{\mathrm{min}}=\frac{3c^2}{a^2+b^2+c^2}1.$$
(7)
For the voids found by the AM algorithm, the denominator of $`Q`$ can be written as

$$\overline{(x^{})^2}+\overline{(y^{})^2}+\overline{(z^{})^2}=\overline{D^2}\overline{D}^2,$$
(8)
where

$$\overline{D}^2=\overline{x}^2+\overline{y}^2+\overline{z}^2$$
(9)
and

$$\overline{D^2}=\overline{x^2}+\overline{y^2}+\overline{z^2}.$$
(10)
The numerator can be written as

$$3\overline{(x^{})^2}=\frac{3}{\overline{D}^2}\left[\overline{D}^4+\overline{x}^2\overline{x^2}+\overline{y}^2\overline{y^2}+\overline{z}^2\overline{z^2}+2\overline{x}\overline{y}\overline{xy}+2\overline{y}\overline{z}\overline{yz}+2\overline{z}\overline{x}\overline{zx}\right].$$
(11)
To compute the distribution of $`Q`$ as a function of void volume $`V`$, we first identify voids using the AM algorithm with a density threshold $`f=0.2`$; to eliminate distortions caused by the artificial boundary conditions at $`z_{\mathrm{max}}=0.12`$, we consider only those voids at redshifts $`z<0.11`$. The mean, $`\mu _Q`$, and the standard deviation, $`\sigma _Q`$, of $`Q`$ for the voids found in this way are plotted in Figure 7, as a function of void volume $`V`$. The values of $`\mu _Q`$ and $`\sigma _Q`$ are computed in bins containing 400 voids apiece. Thus, the expected error in the mean value of $`Q`$ for each bin is $`\sigma _\mu =\sigma _Q/\sqrt{400}=0.05\sigma _Q`$. The value of $`\mu _Q`$ indicates whether voids are preferentially oriented with respect to the line of sight; the value of $`\sigma _Q`$ is a measure of the intrinsic asphericity of the voids.

In the left column of Figure 7, the shape of voids is measured in real space, where there should be no preferential distortions along the line of sight. Indeed, $`\mu _Q`$ is not significantly different from unity for voids of all sizes and at all sampling densities: in all cases, $`|\mu _Q1|<2.5\sigma _\mu `$. There is, however, a significant trend in the standard deviation of $`Q`$. As voids get larger, $`\sigma _Q`$ gets smaller, decreasing from $`\sigma _Q0.34`$ at $`V=5\times 10^8`$ to $`\sigma _Q0.27`$ at $`V>10^6`$. This is a consequence of the fact that large voids are more nearly spherical than small voids.

In the right column of Figure 7, the shape of voids is measured in redshift space, where the distortions due to peculiar velocities may cause systematic distortions along the line of sight. At sufficiently high sampling density (the $`L>L_{}`$ and $`L>L_{}/2`$ samples), large voids are significantly elongated along the line of sight ($`\mu _Q>1`$). In the $`L>L_{}`$ sample, $`\mu _Q>1`$ for volumes $`V>4\times 10^70.2\overline{V}`$. In the $`L>L_{}/2`$ sample, $`\mu _Q>1`$ for volumes $`V>1.3\times 10^70.1\overline{V}`$; the greatest deviation of $`\mu _Q`$ from unity is at $`V=1.25\times 10^6\overline{V}`$, where $`\mu _Q1=6\sigma _\mu `$. We conclude that the distortions along the line of sight caused by peculiar velocities are measurable by the AM algorithm only when the sampling density is sufficiently high (corresponding to a limiting galaxy luminosity of $`L_{}`$ or fainter). For our initial CDM power spectrum, voids with $`V\overline{V}`$ show a significant tendency to be elongated along the line of sight by peculiar velocities.

## 5 Conclusion

In examining the large scale structure of the universe, studies of overdense regions (clusters and superclusters) are usefully complemented by studies of underdense regions, or voids. The properties of voids, such as their volumes and shapes, depend on how voids are defined. Since the literature contains many different void definitions and void-finding algorithms, direct comparison of void properties found in different studies is a risky business.

The VPF, $`P_0(V)`$, depends strongly on the number density $`\nu `$ of galaxies in the survey. Although the UPF, $`P_{0.2}(V)`$, depends less strongly on sampling density than the VPF, it is still true that the characteristic void size $`V_{UPF}`$ depends on the sampling density $`\nu `$, with a higher density yielding a smaller characteristic void size. Thus, in a flux-limited survey, where the mean density of detected galaxies, $`\nu (z)`$, decreases with redshift, using a threshold density $`f\nu (z)`$ will produce a characteristic void size that increases with increasing redshift. This will be an artifact of the decreasing sampling density at large redshift, and does not reflect a change in the underlying large scale structure. In determining how the properties of voids depend on redshift, it is more prudent to extract a volume-limited sample from the data before applying a void statistic like the UPF or a void-finding algorithm.

Of the void-finding algorithms outlined in this paper, the “sphere” algorithm has the virtue of (relative) simplicity. However, its restriction that all voids must be spherical leads to an underestimate of void size and does not permit us to measure the distortions of voids caused by peculiar velocities. Adapting the algorithm of Aikio & Mähönen (1998), we were able to determine more accurately the sizes of underdense voids. As in two-dimensional simulations (Ryden & Melott, 1996), the effect of peculiar velocities is to increase the characteristic void size. The AM algorithm also permits us to measure the elongation of voids along the line of sight. In real space, it is found that large voids are intrinsically more nearly spherical than smaller voids. This can be regarded as a manifestation of the tendency for large voids within a bubbly structure to expand and become more nearly spherical at the expense of their smaller neighbors (Regös & Geller, 1991). In redshift space, large voids are seen, for the CDM spectrum used in our simulations, to have a slight but statistically significant tendency to be elongated along the line of sight. Note, however, that the void distortions can only be detected at a sufficiently high sampling density. Future redshift surveys such as that provided by the Sloan Digital Sky Survey will provide sufficiently high galaxy densities and a large enough number of voids to accurately measure the peculiar velocity distortion of voids in the real universe. When deeper redshift surveys are available, the (relatively small) peculiar velocity distortions can be subtracted out to reveal the cosmological distortions resulting from the deceleration of the Hubble expansion.

This work was supported by a grant from the Ohio Supercomputer Center. ALM wishes to acknowledge support from the National Science Foundation under grant number AST-0070702, the University of Kansas General Research Fund and the National Center for Supercomputing Applications.


format: text

text:
# Behavior of the thermopower in amorphous materials at the metal-insulator transition

## I Introduction

Transport phenomena in disordered quantum systems have been studied for many years,  yet many open problems remain. One focus of these investigations is the metal-insulator transition (MIT). This quantum phase transition from a good conducting material to an insulator may be induced by disorder due to localization  or by interactions such as electron-electron interactions and electron-lattice coupling. In three-dimensional (3D) amorphous materials the MIT is mainly attributed to disorder. For example, in heavily doped semiconductors the disorder is brought about by the random distribution of dopant atoms in the crystalline host. However, indications of electron-electron interactions have also been found, e.g., in the d.c. conductivity $`\sigma `$ (or resistivity $`\rho =1/\sigma `$) in doped semiconductors in both metallic and insulating regimes.

A further open problem is the behavior of the thermoelectric power $`S`$ or the Seebeck coefficient of disordered materials near the MIT. In amorphous alloys and both compensated Si:(P,B) and uncompensated Si:P, $`S`$ continuously changes from negative to positive values or vice versa at low temperature $`T`$. This corresponds to a change of thermal conductors from electrons to holes or conversely and has been attributed to electron-phonon interaction in amorphous alloys.  On the other hand in heavily doped semiconductors the sign change is believed to be caused by electron-electron interactions or attributed to the existence of local magnetic moments and their interactions with electrons.  This conclusion is based on the suppression of the anomalous behavior by a magnetic field.  We remark that the sign change in $`S`$ is also observed in metals, high-$`T_c`$ materials and quasicrystals. Analytical treatments of metals as a degenerate free-electron gas taking into account inelastic scattering with phonons as well as numerical considerations incorporating electronic correlations in superconductors have also been shown to generate a sign change in $`S`$. But in these systems the sign change occurs at a $`T`$ value which is 2 orders of magnitude higher than in disordered systems. Note that $`S`$ is mainly due to two distinct effects: (i) the diffusion of the charge carriers and (ii) the net momentum transfer from phonons to carriers. But for $`T<0.3`$ K as considered in this work, the diffusive part of the thermopower dominates that of the phonon-drag contribution. Hence, from this point on in this paper $`S`$ denotes only the diffusion thermopower.

The prototype for a theoretical description of 3D disordered systems is the Anderson model of localization. Near the MIT at $`T=0`$, $`\sigma `$ behaves as

$$\sigma _\text{c}\{\begin{array}{cc}\sigma _0\left|1\frac{E_\text{F}}{E_\text{c}}\right|^\nu ,\hfill & |E_\text{F}|E_\text{c},\hfill \\ 0,\hfill & |E_\text{F}|>E_\text{c},\hfill \end{array}$$
(1)
where $`E_\text{F}`$ is the Fermi energy, $`E_\text{c}`$ is the mobility edge which separates the extended conducting states from localized insulating states, and $`\nu `$ is a universal critical exponent.  By using Eq. (1) for $`\sigma `$ in a linear response formulation the behavior of the thermoelectric transport properties such as $`S`$, the thermal conductivity $`K`$ and the Lorenz number $`L_0`$  at the MIT have been computed. Moreover, similar to $`\sigma `$, the quantities $`S`$, $`K`$ and $`L_0`$ have also been found to obey scaling. The scaling form of the dynamical conductivity $`\sigma `$ close to the MIT in 3D is given as

$$\frac{\sigma (t,T)}{T^{1/z}}=\left(\frac{t}{T^{\nu z}}\right).$$
(2)
Here $`t`$ measures a dimensionless distance from the critical point, such as $`t=(E_\text{F}E_\text{c})/E_\text{c}`$, the correlation-length exponent $`\nu `$ in 3D is equivalent to the conductivity exponent as given in Eq. (1), and $`z`$ is the dynamical exponent. For a noninteracting system such as the Anderson model, one expects $`z=d`$ in $`d`$ dimensions. But, instead of obtaining $`z=3`$ in the scaling form of $`\sigma `$, one finds $`z\nu =1`$.  In addition to this discrepancy, $`S`$ turns out to be at least one order of magnitude larger  than the experimental results in doped semiconductors  and in amorphous alloys.  Furthermore, the sign change in $`S`$ cannot be explained using the Anderson model and Eq. (1). One may argue that the discrepancies between the transport calculations and the experimental measurements are due to the absence of interactions in the Anderson model. Indeed, interactions may influence the behavior of the thermoelectric transport properties. Yet we emphasize that the neglect of interactions in the Anderson model is not entirely inconsistent with the experimental situation in 3D amorphous materials. For example, recent measurements in Si:P yield $`\sigma `$ scaling with $`z3`$ and $`\nu 1`$. This agrees with $`z=d`$ as predicted by the scaling arguments  for noninteracting systems.

The goal of this paper is to show that the correct value of $`z`$, the right order of magnitude of $`S`$ at the MIT, and perhaps even the sign change in $`S`$ at low $`T`$, can be described within a linear response formulation using the noninteracting Anderson model of localization. However, in order to do so, we have to use a more suitably chosen energy-dependent $`\sigma _\text{c}`$ instead of Eq. (1). After a brief review of linear transport theory, we construct a new form for $`\sigma _\text{c}`$ as a function of energy $`E`$ and $`T`$ from experimental data. By using this model data as input for the linear response formulation, we compute the temperature dependence of $`S`$, $`K`$, $`L_0`$ and also $`\sigma `$ and show that they have the expected qualitative and quantitative behavior close to the MIT. Finally, we show that a small variation in $`\sigma _\text{c}(E,T)`$ can change the sign of $`S`$. This effect cannot be produced simply by varying the density of states $`\varrho `$ or the chemical potential $`\mu (T)`$.

## II Linear Thermoelectric Transport Theory

In the presence of a small temperature gradient $`T`$, the electric current density $`𝐣_1`$ and the heat current density $`𝐣_2`$ induced in a system are given (to linear order) as

$$𝐣_i=|e|^i\left(|e|L_{i1}𝐄L_{i2}T^1T\right),$$
(3)
where $`e`$ is the electron charge and $`𝐄`$ is the induced electric field. $`L_{ij}`$ are the kinetic coefficients. Since we do not consider the presence of a magnetic field in this work, the Onsager relation $`L_{ij}=L_{ji}`$ holds. Ohm’s law,

$$𝐣_1=\sigma 𝐄,$$
(4)
implies that in Eq. (3)

$$\sigma =L_{11}.$$
(5)
The flow of thermal conductors due to $`T`$ is counteracted by an electric force arising from $`𝐄`$ making $`𝐣_1=0`$. Equation (3) then yields the thermoelectric power $`S`$ which relates $`T`$ to $`𝐄`$,

$$S=\frac{L_{12}}{|e|TL_{11}}.$$
(6)
The sign of $`S`$ determines whether the thermal carriers are electrons or holes. Using the Sommerfeld expansion for $`|E_\text{F}E_\text{c}|>k_\text{B}T`$, $`S`$ is given as

$$S=\frac{\pi ^2k_\text{B}^2T}{3|e|}\frac{\text{d}\mathrm{ln}\sigma (E)}{\text{d}E}|_{E=E_\text{F}},$$
(7)
where $`E_\text{F}`$ is the Fermi energy, $`k_\text{B}`$ is Boltzmann’s constant and $`\sigma (E)`$ is assumed to be a slowly varying function on the scale of $`k_\text{B}T`$. Equation (7) is also known as the Mott formula.

The thermal conductivity $`K`$ determines the contribution to $`𝐣_2`$ stemming from $`T`$. Using Eqs. (5) and (6) in $`𝐣_2`$ we obtain $`K`$ in terms of the kinetic coefficient as

$$K=\frac{L_{22}L_{11}L_{21}L_{12}}{|e|^2TL_{11}}.$$
(8)
From the definition of the Lorenz number it follows that

$$L_0\frac{e^2}{k_\text{B}^2}\frac{K}{\sigma T}=\frac{L_{22}L_{11}L_{21}L_{12}}{(k_\text{B}TL_{11})^2}.$$
(9)
In metals at room $`T`$, $`L_0=\pi ^2/3`$.. It also takes on the same value at $`T10`$ K in metals where the electrons suffer no inelastic scattering processes.

The primary consideration then in determining $`\sigma `$, $`S`$, $`K`$ and $`L_0`$ is to calculate $`L_{ij}`$. Under the assumptions that the system is noninteracting and inelastic scattering processes are absent, $`L_{ij}`$ are given in the Chester-Thellung-Kubo-Greenwood formulation  as

$$L_{ij}=_{\mathrm{}}^{\mathrm{}}A(E)[E\mu (T)]^{i+j2}\left[\frac{f(E,\mu ,T)}{E}\right]𝑑E,$$
(10)
where $`i,j=1,2`$, $`\mu (T)`$ is the chemical potential, $`f(E,\mu ,T)`$ is the Fermi distribution function and $`A(E)`$ contains all the system-dependent features.

Lastly, we note that the $`T`$ dependence of $`\mu `$ can be obtained for noninteracting systems from

$$n(\mu ,T)=_{\mathrm{}}^{\mathrm{}}𝑑E\varrho (E)f(E,\mu ,T)$$
(11)
where $`n`$ is the number density of electrons and $`\varrho `$ is the density of states. Knowing $`\varrho `$ and keeping $`n`$ constant, we find numerically that $`\mu (T)T^2`$ in the 3D Anderson model with an increased effective mass due to the disorder.

## III A Phenomenological Approach

There are only two parameters that are model dependent in the transport theory discussed in Sec. II. These are $`A(E)`$ and $`\mu (T)`$. In order to determine the behavior of the thermoelectric transport properties close to the Anderson MIT, $`A(E)`$ in Eq. (10) has usually been set  equal to the critical behavior of $`\sigma `$ given by Eq. (1). As mentioned in the introduction, this leads to the unphysical value for $`z=1/\nu `$ and therefore an unphysical frequency and $`T`$ dependence of $`\sigma `$. The main reason for this behavior is easily understood: there is no $`T`$ dependence in Eq. (1) and thus all $`T`$ dependence in Eq. (5) is due to the broadening of the Fermi function in Eq. (10) with increasing $`T`$. Thus in order to model the correct $`T`$ dependence, we should add to Eq. (1), valid at $`T=0`$, the desired $`T`$ dependencies such as $`\sigma T^{1/z}`$ in the metallic and $`\sigma \mathrm{exp}(T)`$ in the insulating (say, variable-range-hopping) regimes. Such a purely theoretical model for $`\sigma _\text{c}(E,T)`$ will then incorporate a multitude of constants that can be adjusted to fit the experimental results. Of course this is of limited practical use since the validity of these fitting parameters is hard to justify.

Here we will instead use as input for $`\sigma _\text{c}(E,T)`$ recent experimental data obtained by Waffenschmidt et al.  who measured $`\sigma `$ in Si:P at the MIT under uniaxial stress. Their data yield good scaling of $`\sigma `$ according to Eq. (2) with a dynamical exponent $`z=2.94\pm 0.3`$ and $`\nu =1\pm 0.1`$. These values agree with the scaling arguments  and reasonably well with the numerical results  for noninteracting systems.

In the $`\sigma (t,T)`$ scaling of Ref. , $`t=(ss_\text{c})/s_\text{c}`$ where $`s`$ is the stress and $`s_\text{c}`$ the corresponding value at the transition. We sample their scaled data for several values of $`(t,T)`$ and fit a spline curve  $`\sigma _\text{c}`$ to these points in order to get a smooth functional form for $`\sigma (t,T)`$. Transforming the spline $`\sigma _\text{c}`$ as a function not only of $`T`$ but also of $`E`$, we set $`t=(EE_\text{c})/E_\text{c}`$. Finally, we substitute $`\sigma _\text{c}(E,T)`$ for $`A(E)`$ in Eq. (10). Consequently, the thermoelectric transport properties are directly obtained from Eqs. (5)–(9).

In this paper we consider temperatures from $`0.01`$ K to $`0.2`$ K. Far from the transition we could not probe lower than $`T<0.02`$ K. This is due to the limited input data and consequently a limited range of the spline function that generated $`\sigma _\text{c}(E,T)`$. The unit of $`\sigma _\text{c}`$ is taken as $`\mathrm{\Omega }^1\text{cm}^1`$ consistent with the experiments. The $`E`$ scale is (arbitrarily) fixed at $`1`$ meV which is the order of magnitude of the binding energy of an isolated donor in a heavily doped semiconductor.  In order to compare with the previous results in the Anderson model  we let $`E_\text{c}=7.5`$. We emphasize that this value is of no significance and can be assigned (nearly) arbitrarily. The important point to consider is the location of the Fermi energy $`E_\text{F}`$ with respect to $`E_\text{c}`$. This distinguishes the electronic regimes. Thus, the metallic, critical and insulating regimes are identified as $`|E_\text{F}|<E_\text{c}`$, $`E_\text{F}=E_\text{c}`$, and $`|E_\text{F}|>E_\text{c}`$, respectively. Usually, $`\mu (T)`$ is derived in Eq. (11) from $`\varrho `$ of the 3D Anderson model of localization. In the next section we shall also show the effect of using a different functional form of $`\mu (T)`$.

## IV Results and discussions

### A Temperature dependence of the thermoelectric transport properties

Consistent with the dynamics of the experiment in Ref. , we expect $`\sigma (T)T^{1/z}`$ at the critical regime with $`z3`$. This is indeed the behavior of $`\sigma (T)`$ close to $`E_\text{c}`$ as we show in Fig. 1. For $`|E_\text{F}E_\text{c}|0.2`$ meV we obtain $`z=3.2\pm 0.3`$. Note that $`\sigma (T)=L_{11}`$ has been integrated according to Eq. (10) over the energy range where $`f/E10^{20}`$ meV. Thus our numerical calculation of $`\sigma `$ is consistent since it reproduces closely the original result in Ref. . If we plot the results in Fig. 1 with respect to $`(\mu E_\text{c})/E_\text{c}T^{1/\nu z}`$ we obtain a rough scaling of $`\sigma `$ similar to Fig. 3 of Ref. .

We next turn our attention to the thermoelectric power $`S`$. In the 3D Anderson model of localization, we know that when using Eq. (1) one obtains $`S0`$ in the metallic regime  while in the insulating regime $`S`$ does not approach zero but seems to diverge as $`T0`$.  At the MIT $`S`$ is a constant  of the order of $`100\mu `$V/K. In Fig. 2, we show that in the present approach $`S`$ in the vicinity of the MIT is two orders of magnitude smaller compared to previous results for the Anderson model. The magnitude of $`S`$ is in fact comparable to the experimental results in disordered systems. Furthermore, $`S0`$ as $`T0`$ in the metallic, critical and insulating cases. This behavior of $`S(T)`$ in all electronic regimes was observed  in amorphous AuxSb1-x and in amorphous Ge1-xAux. As indicated by the different lines in Fig. 2, $`S`$ is in good agreement with Eq. (7) since $`\sigma `$ in Fig. 1 is smooth across the transition at finite $`T`$. Note that in order to evaluate Eq. (7) properly for the system considered here, the $`E`$ dependence of the input spline $`\sigma _\text{c}(E,T)`$ was used instead of $`\sigma (T)`$ from Fig. 1. We emphasize that it is no contradiction that $`S`$ is positive here but mainly negative in the doped semiconductors in all electronic regimes. In the energy regions close to $`E_\text{c}>0`$ the charge carriers are holes instead of electrons as shown in Ref. . $`S`$ would be negative if we had chosen the left mobility edge $`E_\text{c}<0`$ for low filling.

The corresponding $`T`$ dependence of $`K`$ is shown in Fig. 3. We find that $`K0`$ as $`T0`$ in all electronic regimes. This is also the behavior of $`K`$ using $`\sigma _\text{c}`$ in Eq. (1). In the metallic regime $`K`$ is larger than in the insulating regime since there are more heat carriers in the former. From the results of $`\sigma `$ and $`K`$ in Figs. 1 and 3 respectively, we obtain $`L_0`$. As shown in Fig. 4, $`L_0\pi ^2/3`$ as $`T0`$ whether it be in the metallic, critical or insulating regime. This is different from the results using Eq. (1) for $`\sigma `$. There one obtains an $`L_0`$ which depends on the conductivity exponent in the critical and insulating regimes while it approaches the universal value  $`\pi ^2/3`$ only in the metallic regime. Here we see no markedly distinct behavior in the metallic regime compared to the insulating regime. For $`|E_\text{F}E_\text{c}|=0.1,0.2`$ meV, $`L_0`$ in the metallic regime is less than its corresponding value in the insulating regime. For $`|E_\text{F}E_\text{c}|=0.5,1.0`$ meV, $`L_0`$ in the metallic regime is larger than its corresponding value in the insulating regime.

In the calculation of $`\sigma (T)`$, $`S(T)`$, $`K(T)`$ and $`L_0(T)`$, we used a phenomenological construction of $`\sigma _\text{c}(E,T)`$. Furthermore, we have assumed that the density of states $`\varrho `$ is the same as that of the 3D Anderson model of localization given in Ref. . Since this $`\varrho (E)`$ is a smooth and (restricted to $`E>0`$) monotonic function, $`\mu (T)`$ obtained from Eq. (11) is also smoothly and monotonically varying with $`T`$ as described in Sec. II.

### B Effects of a structured $`\varrho `$

We now consider the effects of a possible structure in $`\varrho `$. We shall assume here that this structure corresponds only to variations in $`\mu (T)`$ and not in $`\sigma `$. In Fig. 5 we show two examples of a modified $`\mu (T)`$ in the critical regime. Example A has a pronounced maximum, while example B has both a maximum and a minimum. The height of the maximum in both examples A and B is $`0.1`$ meV. Note that this is significantly larger than the halfwidth of the bump which is $`0.005`$ meV. This is also true for the depth of the minimum in example B. Thus a small change in $`T`$ corresponds to a large change in $`\mu (T)`$. Applying these forms of $`\mu (T)`$ together with Eq. (1) for $`\sigma `$ reproduces the same structures in $`S`$. For example, using form B of $`\mu (T)`$ we obtain an $`S`$ having both a maximum and a minimum in the same $`T`$ interval as $`\mu (T)`$. But $`S`$ is still of the order of $`100\mu `$V/K while the variations are only of the order of $`10\mu `$V/K and not large enough to cause a sign change in $`S`$. On the other hand using the phenomenological construction of $`\sigma _\text{c}(E,T)`$ yields even smaller changes. We observe variations in $`L_{11}`$ and $`L_{12}`$ of less than 10$`\%`$ from their unmodified values. Consequently, we find negligible changes in $`S`$. Figure 2 would appear unmodified. Hence we conclude that even a large change in the density of states $`\varrho `$ and thus also in $`\mu (T)`$ is not sufficient to cause the change of sign for $`S`$ as observed in the experiments. However, this weak dependence of $`S`$ on $`\varrho `$ and $`\mu (T)`$ at least justifies our use of the simple Anderson density of states in the present paper.

### C Effects of a structured $`\sigma (E)`$

Let us now assume that for small $`T`$ there are non-monotonicities in $`\sigma _\text{c}(E,T)`$ — although these have not been observed in the experiments. Thus we consider the case when there is a sizeable change in $`\sigma _\text{c}(E,T)`$ in the region close to $`E_\text{c}`$ for small $`T`$. The corresponding ‘bumps’ in $`\sigma _\text{c}(E,T)`$ are shown in Fig. 6 with different peak heights and with half-widths $`<1`$ meV. For simplicity they are essentially quadratic functions of $`E`$ and have been generated such that they decay quickly as $`\mathrm{exp}(T^4)`$ with increasing $`T`$. The height of the bumps is $`<1\mathrm{\Omega }^1`$cm-1 which is at least an order of magnitude smaller than the values observed for $`\sigma `$ in the measurements. The lowest temperature studied is $`T=6`$ mK and we shall only consider metallic and insulating regions with $`|E_\text{F}E_\text{c}|0.1`$ meV.

Our results in Fig. 7 using $`\sigma _\text{c}(E,T)`$ with and without bumps indicate that there are only small variations in the slope of $`\mathrm{log}(\sigma )`$ and $`\sigma T^{1/3}`$ remains valid within the accuracy of these estimations. We remark that the lowest measured temperature in Ref.  is $`15`$ mK. From Fig. 7 we see that the variations for $`T15`$ mK are much smaller than for $`T<15`$ mK. Hence, these variations in $`\sigma `$ could not have been observed in the experiments.

In Fig. 8 we show how the bumps affect $`S`$. Even with the very small bump 3, $`S`$ changes sign in the critical regime as $`T0`$. As the bump increases this change becomes more pronounced. The temperature $`T_{S=0}`$ at which the sign change occurs is $`0.1`$ K for bump 3, $`T_{S=0}=0.2`$ K for bump 2, and $`T_{S=0}0.4`$ K for bump 1. These results for $`T_{S=0}`$ are still about one order of magnitude less than in semiconductors  and two orders of magnitude smaller than in amorphous alloys. Of course, as shown in Fig. 8, $`T_{S=0}`$ shifts to higher values as the bump height increases. Nevertheless, the minimum value of $`S`$ for $`T<T_{S=0}`$ has the same order of magnitude as the corresponding maximum value of $`S`$ in Si:(P,B) and in Si:P. We emphasize that the value of $`T_{S=0}`$ of course depends on the energy unit chosen and thus will vary for systems with different bandwidths, e.g., a larger band width will give rise to a larger value of $`T_{S=0}`$. We observe a similar sign change in the metallic regime but the depth of the minimum is smaller than in the critical regime. The Mott formula (7) with $`\sigma (T)`$ given in Fig. 7 can readily model this behavior since $`\sigma (T)`$ remains slowly varying even if $`\sigma _\text{c}(E,T)`$ has a bump near $`E_\text{c}`$. In the insulating regime, $`S`$ has a shallow maximum and drops back to zero as $`T0`$. This is different from experiment, where $`S`$ changes sign and neither has a maximum nor minimum in the insulating regime.

Unlike in $`S`$, there is no dramatic change in $`K`$ as can be seen in Fig. 9. We find only negligible variations at $`T<15`$ mK. This should be expected since there has also been hardly any change in the slope of $`\sigma (T)`$ except for $`T<15`$ mK. However, the small increase in $`\sigma `$ at $`T<15`$ mK in the metallic regime together with the minimally modified $`K`$ leads to a drastic change in $`L_0`$ even in the case for the smallest bump. The increase and decrease in $`\sigma `$ leads to a maximum and minimum in $`L_0`$, respectively. However, $`L_0`$ still approaches the universal value $`\pi ^2/3`$ for $`T0`$ as demonstrated in Fig. 10.

## V Conclusions

In this paper, we have shown that the anticipated value of the dynamical scaling exponent $`z3`$ as well as the right order of magnitude for the thermopower $`S1\mu `$V/K at the MIT can be obtained when taking into account the expected $`T`$ dependencies in addition to the simple scaling behavior of Eq. (2). Our approach is phenomenological in the sense that we have refrained from using fitting parameters and have rather taken experimental data as input. Using this data, we can explain the large deviations from experimental results as reported in the theoretical studies of Refs. . We have shown that our results for $`S`$ agree with those predicted by the Mott formula (7) since we have used a $`\sigma `$ slowly varying on the scale of $`k_\text{B}T`$ near the MIT. We emphasize, however, that for a disordered system where interactions are negligible, we should still expect the Anderson-type transition as given in Eq. (1) at $`T=0`$. Consequently, $`S100\mu `$V/K at the MIT  should again be expected and one should observe a large increase of $`|S|`$ at very low $`T`$. However, such temperatures appear presently inaccessible by experiment.

As a further challenge, we considered the sign changes observed in $`S`$ at low $`T`$. We found that even large variations in the chemical potential $`\mu (T)`$ do not lead to a sign change in $`S`$. On the other hand, a variation in the input $`\sigma _\text{c}(E,T)`$ data can give rise to such a sign change in $`S`$, while at the same time resulting in only small changes in the conductivity $`\sigma `$. Hence we have effectively modeled the underlying physical reasons for the sign change — which have been attributed to electron-electron interactions or to the existence of local magnetic moments and their interactions with electrons or to inelastic scattering with phonons — by simply changing the input $`\sigma _\text{c}(E,T)`$. Regarding a possible test for the existence of such a structured $`\sigma _\text{c}(E,T)`$, we have shown that the $`T`$ variation of $`L_0`$ is much more sensitive to the bumps than $`\sigma `$. Thus we have been able to describe the main features of the critical behavior of $`S(T)`$ although it remains unclear what might cause bumps in $`\sigma _\text{c}(E,T)`$ close to $`E_\text{c}`$. A microscopic and possible system-dependent approach to the problem may eventually account for these abrupt changes in $`\sigma _\text{c}`$. Of course, if many-particle interactions and electron-phonon coupling are important, we no longer expect the feasibility of the Chester-Thellung-Kubo-Greenwood formulation  used here.

## Acknowledgment

This work was supported by the DFG through Sonderforschungsbereich 393, by the DAAD, the British Council and the SMWK. The authors would like to thank R. Fletcher, R. Rentsch, R. Rosenbaum and B. Sandow for useful discussions. We thank R. Rosenbaum for a critical reading of the manuscript.