text:
# Elastic properties of the vortex lattice for a superconducting film of finite thickness

## Abstract

In this paper we investigate the elastic properties of the vortex lattice for a superconducting film of finite thickness. We derive an analytic expression for the compression modulus. The shear modulus is evaluated numerically by using both the Pearl interaction potential, valid in the limit of very thin film, and a potential for films of arbitrary thickness. A comparative study of the shear moduli is carried out.

The problem of a vortex emerging perpendicularly to a surface of a superconductor was first considered by Pearl. He pointed out that the vortex-vortex interaction potential at the surface of superconductors follows the power law $`1/r`$ ($`1/k`$ in Fourier space) at large distances. The Pearl interaction potential has been widely used to investigate the equilibrium and the elastic properties of the vortex lattice in superconducting films.

In this work we will show that in fact the Pearl interaction potential is capable of describing satisfactorily the superconducting properties of a very thin superconductor. However, we will show that as the thickness of the film approaches the London penetration depth $`\lambda `$, the use of this potential may underestimate significantly the value of the shear modulus for sufficiently low induction.

It has been shown elsewhere that the energy of an ensemble of interacting vortices inside a superconducting film of finite thickness $`d`$ and the energy of the stray fields in the vacuum is given by

$$F=\frac{\mathrm{\Phi }_0^2}{8\pi }\frac{d^2k}{(2\pi )^2}\frac{1}{\lambda ^2\alpha ^2}\left\{d+\frac{2}{\lambda ^2k\alpha \left[\alpha +k\mathrm{coth}\left(\frac{\alpha d}{2\lambda }\right)\right]}\right\}|S(𝐤)|^2.$$
(1)
where $`𝐤=(k_x,k_y)`$, $`k=\sqrt{k_x^2+k_y^2}`$, $`\alpha =\sqrt{k^2+\lambda ^2}`$, and $`\mathrm{\Phi }_0`$ is the quantum flux. The structure factor is given by

$$S(𝐤)=\underset{i}{}e^{i𝐤𝐑_i}.$$
(2)
Here $`𝐑_i`$ is the position of the $`i`$th vortex line. Note that this result should be valid for an ensemble of distorted vortices, that is, the positions of the vortices do not necessarily correspond to the equilibrium configuration. The first term inside Eq. (1) represents the interaction energy of the vortex lines as if the surfaces were absent. The second term represents the surface energy associated to the magnetic energy of the stray field at the vacuum. Note that for $`k`$ small (large $`r`$), $`\alpha ^21/\lambda ^2`$. Thus, the surface energy goes as $`\mathrm{\Phi }_0^2/8\pi ^2r`$. This is the Pearl result for vortices emerging from a semi-infinite isotropic superconductor. Another interesting particular case of Eq. (1) is the limit of a very thin film, $`d0`$, and $`k`$ small. In this limit, from Eq. (1) it is straightforward to show that

$$F=E_0\frac{dk^2}{(2\pi )^2}\frac{2\pi d}{k\mathrm{\Lambda }^1+k^2}|S(𝐤)|^2,$$
(3)
where $`E_0=(\mathrm{\Phi }_0/4\pi \lambda )^2`$, $`\mathrm{\Lambda }=2\lambda ^2/d`$. This is precisely the energy of an ensemble of interacting vortices in a very thin film first obtained by Pearl. However, his derivation is supposed to be valid for any $`k`$. In this work we use Eq. (3) with no restriction on $`k`$.

The integrand of Eq. (1) contains two terms. The bulk term follows the power law $`1/k^2`$ for large $`k`$. The integrand of Eq. (3) has the same behavior as its corresponding one. As a result, the self-energy contributions to both Eqs. (1) and (3) diverge. The divergence is logarithmic. This is so because London theory neglects the size of the vortex-core. To remove this short-length scale divergence we use the Gaussian cutoff. This regularization procedure consists in multiplying the integrand of Eq. (1) and (3) by a factor $`e^{2\xi ^2k^2}`$, where $`\xi `$ is the coherence length. In the first case, only the bulk term needs a cutoff. The London picture is valid in the limit of low induction and very strong type II superconductor, $`\kappa =\lambda /\xi 1`$.

Let us suppose that initially the vortex lines are located at the equilibrium positions $`𝐑_i^0=(X_i^0,Y_i^0)`$. Let us denote by $`𝐮_i𝐮(𝐑_i^0)`$ the displacement of the $`i`$th vortex line from its equilibrium position. The new position of the $`i`$th vortex line is then given by $`𝐑_i=𝐑_i^0+𝐮_i`$. The elastic energy of the distorted vortex lattice is most conveniently written in $`𝐤`$-space by expressing the displacement in terms of its Fourier transform

$$𝐮_i=_{BZ}\frac{dk^2}{(2\pi )^2}\frac{1}{d}𝐮(𝐤)e^{𝐤𝐑_i^0},$$
(4)
where the integration is taken over the first Brillouin zone (BZ) of the vortex lattice. Upon taking this into Eq. (1) and expanding it up to second order in $`𝐮_i`$ we obtain for the excess free energy

$$\delta F=\frac{1}{2}\frac{dk^2}{(2\pi )^2}\frac{1}{d}u_\alpha (𝐤)\mathrm{\Phi }_{\alpha \beta }(𝐤)u_\beta (𝐤),$$
(5)
with $`(\alpha ,\beta )=(x,y)`$. Here $`B=\mathrm{\Phi }_0/A`$ is the induction, $`A`$ is the area of a unit cell of the vortex lattice. The factor $`1/d=_{\pi /d}^{\pi /d}𝑑k_z/(2\pi )`$ in the integrals was introduced deliberately to resemble the three dimensional corresponding expressions. The coefficients $`\mathrm{\Phi }_{\alpha \beta }(𝐤)`$, called elasticity matrix, are real, symmetric, and periodic in $`𝐤`$-space. One has

$$\mathrm{\Phi }_{\alpha \beta }(𝐤)=\frac{B^2}{4\pi }\underset{𝐐}{}[(𝐤+𝐐)_\alpha (𝐤+𝐐)_\beta V(𝐤+𝐐)Q_\alpha Q_\beta V(𝐐)],$$
(6)
where $`𝐐`$ are the reciprocal lattice vectors. Although the last equation is valid for any type of vortex lattice, we will use a triangular lattice. The basis vectors of the reciprocal lattice are given by

$$𝐐_1=\frac{2\pi }{a}\left(\widehat{𝐱}\frac{\widehat{𝐲}}{\sqrt{3}}\right),𝐐_2=\frac{4\pi }{a}\frac{\widehat{𝐲}}{\sqrt{3}}$$
(7)
where $`a^2=2\mathrm{\Phi }_0/\sqrt{3}B`$. Then, the vortex positions in reciprocal space are given by $`𝐐𝐐_{mn}=n𝐐_1+m𝐐_2`$, with $`m,n`$ integers.

The interaction potential $`V(𝐤)`$ for a film of arbitrary thickness is

$$V(𝐤)=\frac{1}{\lambda ^2\alpha ^2}\left[1+\frac{2}{d\lambda ^2k\alpha (k+\alpha )}\right],$$
(8)
and in the limit of Pearl is

$$V(𝐤)=\frac{2}{kd+\lambda ^2k^2}.$$
(9)
Since we have assumed that the vortex lines are straight and parallel to each other, there will be only two elastic constants connected to the excitation of the lattice. One is the compression modulus and the other one is the shear modulus. Usually, shear is softer than compression. Let us first neglect the discreteness of the lattice by considering only the $`𝐐=0`$ contribution to Eq. (6). Within this continuum approximation, only the compression modulus is obtained. One has

$$c_{11}(𝐤)=\frac{B^2}{4\pi }V(𝐤).$$
(10)
Note that in the long wave length limit, $`k0`$, the interaction potential diverges, both for the full expression and within the Pearl limit \[see Eqs. (8) and (9)\]. As a result, in this local limit the compression modulus diverges. Thus, as the vortex density becomes very small, the energy cost to compress the lattice is very large.

The evaluation of the more important shear modulus is possible only if one goes beyond the continuum limit. Over the first Brillouin zone the shear modulus is nearly a constant, that is, it does not contain any significant non-locality like the compression modulus. Thus, either of the following definitions is sufficient to describe the shear deformations

$$c_{66}=\underset{k0}{lim}\frac{\mathrm{\Phi }_{yy}(k,0)}{k^2},c_{66}=\underset{k0}{lim}\frac{\mathrm{\Phi }_{xx}(0,k)}{k^2}.$$
(11)
In what follows, the induction $`B`$ is in units of the upper critical field $`H_{c2}`$, $`b=B/H_{c2}`$; lengths are in units of $`\lambda `$; the Ginzburg-Landau parameter used is $`\kappa =50`$. By using these parameters we evaluated numerically $`c_{66}^p`$ and $`c_{66}`$, where the superscript $`p`$ means that the Pearl potential was used to calculate the shear modulus, whereas with no superscript we used the full expression of Eq. (8). Fig. 1 shows the ratio $`c_{66}^p/c_{66}`$ as a function of the film thickness. As one can see from this figure, for sufficiently low induction the difference between both shear modulus grows as $`d`$ increases. Note that this occurs for values slightly larger than $`\xi `$. However, these differences tend to disappear as the induction increases.

In Fig. 2 we also show a plot of the ratio $`c_{66}^b/c_{66}`$ as a function of the film thickness, where the superscript $`b`$ stands for bulk. Note that as the induction increases the shear modulus the bulk shear modulus becomes approximately equal to full shear modulus. This was first experimentally observed by Fiory.

The conclusion we draw from this scenario is that the Pearl potential is indeed valid in the limit of very thin film. However, the study of the elastic properties of the vortex lattice in superconducting films, of thickness slightly larger than $`\xi `$ and the induction very close to the lower critical field $`H_{c1}`$, might become non-reliable by employing this potential. In this situation it should be more convenient to use the full expression for the interaction potential, but for films of thickness not much larger than $`\lambda `$, because in this limit bending of the vortex lines may become important.

###### Acknowledgements.

The author thanks the Brazilian Agencies FAPESP and CNPq for financial support.


format: text

text:
# Anomalous Diffusion in Quasi One Dimensional Systems

## Abstract

In order to perform quantum Hamiltonian dynamics minimizing localization effects, we introduce a quasi-one dimensional tight-binding model whose mean free path is smaller than the size of the sample. This size, in turn, is smaller than the localization length. We study the return probability to the starting layer using direct diagonalization of the Hamiltonian. We create a one dimensional excitation and observe sub-diffusive behavior for times larger than the Debye time but shorter than the Heisenberg time. The exponent corresponds to the fractal dimension $`d^{}0.72`$ which is compared to that calculated from the eigenstates by means of the inverse participation number. PACS numbers: 72.10.Bg, 73.20.Dx, 73.20.Fz

,

The classical kinetic theory predicts that, in a disordered system, the return probability of an excitation decays with a diffusive law $`P(t)(4\pi Dt)^{d/2}`$, where $`D`$ is the diffusion coefficient and $`d`$ the dimension of the system. This fact has been extensively used in many areas of physics, in particular in electronic transport. However, the theory of quantum localization, on the basis of steady state transport properties, made clear that different regimes arise as a function of the disorder $`W`$. For $`W\mathrm{}`$, the eigenstates are completely localized and wave packets do not move through the system (insulating phase). If $`W0`$, the motion is completely ballistic with a velocity $`Ja/\mathrm{}`$. Between these two limits the diffusive behavior (metallic phase) is observed. The particles move around freely between collisions for a certain average length $`\mathrm{}`$, called the mean free path.

More recently, theoretical studies have considered the regime where the critical amount of disorder $`W_C`$ is such that the system is at the metal-insulator transition ($`W_C`$ $`zJ`$ is the typical exchange energies with $`z`$ neighbors at distance $`a`$). In that case, the dynamics of the system is diffusive, but with a smaller exponent implying $`d^{}<d`$ . This reduction in the effective dimension is attributed to the (multi)fractality of the eigenstates at the transition. On the other hand, the diffusion of a spin excitation was directly observed  in NMR experiments. When the spin network has a cubic structure filling the space, the intensity of the excitation decays as $`(tt_0)^{1.5}`$, which is the expected value for diffusion in a three dimensional system. However, the same experiment performed in a chain-like structure (powdered polyethylene) shows anomalous exponents in the diffusion of the excitation; namely $`0.9`$ and $`0.7`$ for the crystalline and amorphous parts of the sample respectively. In these systems, the role of disorder is played by typical energy differences $`W`$ between states, which are smaller than exchange energies $`zJ.`$ This assures a diffusive dynamics with $`DJa^2/\mathrm{}`$. Therefore, the effective dimensions $`d^{}=1.8`$ and $`1.4`$ corresponding to these values reflect the spatial structure of the spin network.

In this work we study the conditions to reach the diffusive regime from actual quantum dynamics, namely the numerical study of model Hamiltonians. In particular, we developed a quasi one dimensional tight binding model system, which we called the Stars necklace model, whose basic unit (layer) is a highly connected cluster (see inset in figure 1) with $`N`$ sites and intralayer hopping $`J=V/N^{0.5}`$. Disorder is introduced through on site energies characterized by a random distribution of width $`W`$ . The initial wave function is a packet defined in one layer of the system. We are interested in the probability of return to the layer, $`P`$, which is the sum of probabilities of finding the particle in every site of the initial layer. To calculate the dynamics of the system we perform an exact diagonalization of the Hamiltonian. One key aspect of the numerical calculation is that the mean free path $`\mathrm{}`$, the size of the system $`L`$ and the localization length $`\xi `$ must obey some restrictive relationships. The condition for a diffusive regime is $`\mathrm{}`$ much smaller than $`L`$, hence assuring that the particle will collide many times before reaching the boundaries of the system. In turn, to stay away from the localized regime, $`L`$ must be smaller than $`\xi `$. For a strictly one dimensional wire $`\xi =2\mathrm{}`$, while for strips and bars with a given number of transverse modes (channels) $`M`$ the localization length is expected to go as $`\xi =2M\mathrm{}`$. In our model the $`M=N1`$ channels available for transport have the same group velocity $`v=2Va/\mathrm{}`$. This striking feature of the model allows to reach a one dimensional diffusive dynamics when $`N\mathrm{}.`$ For finite $`N`$ it provides an optimal representation for one-dimensional excitations. We studied many system sizes and amounts of disorder, a typical evolution is shown in figure 1. The particular system for this figure has $`N=12`$ and a perimeter of $`100a`$, $`W=3V`$, $`\mathrm{}6a`$ and $`\xi 140a`$ . We see that after a ballistic time $`\mathrm{}/v`$, the evolution follows a power law which indicates a diffusive behavior. Nevertheless, the exponent of this power law is somewhat different from the expected one dimensional value $`0.5`$. As in the examples mentioned above, this anomalous exponent could be due to a fractal effective dimension of the system. The fitting of the evolution to a power law with a free exponent resulted in an effective dimension of the system $`d^{}0.7`$. In our model, a possible cause for a fractality in the eigenfunctions of the system is disorder. Strongly localized stated in the band tails are confined around some random points. This means that they represent “holes” in the real space allowed to extended wave functions, thus making the effective dimension of the system smaller than the real one. For times longer than the ones showed in the figure, the autocorrelation function saturates, this is a finite size effect (the saturation value depends linearly on the system size). We also observed that a magnetic field does not change the exponent in the power law noticeable, but reduces the value of the saturation, meaning that there are fewer localized states. Another way to study how the eigenstates of energy $`\epsilon `$ are occupying some fraction of the volume of space is by means of the inverse participation number $`p^1=_i\left|\phi _i(\epsilon )\right|^4`$ . For plane waves one obtains that $`p=L^d`$, i.e. equals the volume of the system. For a localized state $`p`$ is proportional to the volume in which the state has a non-vanishing amplitude. However if the states are extended but fractal, in the thermodynamic limit it diverges as $`p=L^d^{}`$, with $`d^{}`$ an effective dimension that may be different from the dimension $`d`$ of the ordered system. We calculated the inverse participation number for each of the eigenstates of the system and through it the effective dimension ratio $`d^{}/d`$ of everyone of them. The results (shown in figure 2) are in very good agreement with the effective dimension calculated through the fitting of the autocorrelation function depicted in Fig. 1. Summarizing, we have introduced a numerical Hamiltonian model whose exact solution shows a regime with sub-diffusive behavior. Moreover, we presented hints of a fractal dimension of the extended eigenstates induced by the presence of disorder. By hindering particles from a fraction of the available real space, disorder induces a weak breaking of the ergodicity that anticipates the non-ergodicity associated with full localization.


text:
# Quantum Field Symbolic Analog Computation: Relativity Model

## 1 Introduction

Recently, there has been considerable interest in what is called quantum computation . The efforts in this regard are to seek improved ways of performing computations or building new types of computing machines. It is hoped that not only faster computation will be achieved , but that better understanding of computational complexity will come about .

Quantum computers might appear to be discrete systems, such as a finite collection of spins or qubits . But this is not the only possibility. Instead of thinking of a computer based on a physical system which is understood in terms of non-relativistic quantum mechanics with spin added on 111 Spin used to be added in a ad hoc fashion to non-relativistic quantum mechanics until Dirac produced his relativistic equation for the electron ., one can also consider a system based on relativistic quantum mechanics .

In another direction, the classical discrete digital computer  has been generalized to include the possibility of computing over the continuum . This is because the latter way of computing over the continuum is more appropriate to the way we do analysis, physics, and engineering problems. So this is a computing model which is based on classical mechanics. But classical mechanics could also be extended to include relativity, getting relativistic mechanics .

Because, in studying atomic phenomena, classical mechanics has been replaced by quantum mechanics, we could also think of more general models of computing  based on adding relativity to quantum theory to get relativistic quantum field theory.

It is natural to consider our physical or quantum systems in the continuum limit. In fact Isaac Newton , when studying gravitation, found it natural to consider a continuous distribution of matter to model the earth’s gravitational action at external points. From continuum quantum mechanics, by combining relativity, we have quantum field theory.

## 2 Generalized Quantum Computation

The continuum limit (in space-time configuration) of quantum mechanics together with the incorporation of Einstein’s special relativity leads naturally to the relativistic quantum theory of fields. By taking the limit as the velocity of light $`c\mathrm{}`$ we expect to get non-relativistic quantum mechanics. The limit as Planck’s constant $`h0`$ gives classical mechanics . Quantum field theory includes all of quantum mechanics, classical mechanics, and much more 222Examples are, discrete anti-unitary symmetries, CPT invariance and the spin statistics connection .. Included also will be unitary transformations and superposition of amplitudes, which are regarded as prerequisites for quantum computation.

It has been possible to generalize quantum computation to relativistic quantum field computation in a certain model . We expect other forms of quantum computation, namely, those based on non-relativistic quantum mechanics, topological quantum field theories  and classical mechanics  to be related to this generalization. The relationship should shed light, not only on computing possibilities , but also on quantum field theory itself.

## 3 Quantum Fields

Non-relativistic quantum mechanics is not complete, because radiative corrections have to be made to it, using field theory. In dealing with a system corresponding to an infinite number of degrees of freedom, it is well known historically that formulations of quantum field theory like perturbation theory lead to infinities resulting in the need for renormalization. Nevertheless, quantum electrodynamics has turned out to be, “the most accurate theory known to man” 333This statement is attributed to Feynman.. Dirac, Schwinger and Feynman are some of the principal contributers to quantum electrodynamics 444The spectacular history of this is related in  and hence to quantum field theory . Relativistic covariance is of paramount importance in correctly performing the renormalization process.

It is useful here to work within the Wightman formulation  of quantum field theory 555The fruitfulness and utility of this formulation, from a current perspective, is available in . We are dealing with fields in the Heisenberg picture, without using perturbation theory, nor any particular time frame related Hamiltonians. The theory is in terms of analytic functions (Wightman functions) of several complex variables. These functions arise from their boundary values which are vacuum expectation values of the form

$$𝒲_m(x_1,x_2,\mathrm{}x_m)=(\mathrm{\Omega },\varphi _1(x_1)\varphi _2(x_2),\mathrm{}\varphi _m(x_m)\mathrm{\Omega })$$

of products of $`m`$ quantum field operators in a separable Hilbert space. The field operators transform according to appropriate unitary spin representations of the Poincaré (inhomogenous $`SL(2,)`$ group. Quantum fields are uniquely reconstructed from these analytic functions by Wightman.

Let the ($`m`$-point) Wightman functions be denoted by $`W(n;z)`$ where $`z`$ denotes the set of $`n`$ complex variables. Here, $`n=sm`$ where $`s2`$ is the space-time dimension; space-time will consist of 1-time and $`(s1)`$-space dimensions 666Because it is not known, at present, how to physically understand concepts like closed time-like loops in more than one time dimension, there will be only one dimension in time .

Because these analytic functions are fundamental to the theory, one is led to computation of holomorphy domains for these functions over the space of several complex variables, $`^n`$ . The mass spectrum is assumed to be reasonable, in the sense that momentum vectors $`p^\mu `$ lie in the closed forward light cone, with time component $`p^0>0`$ except for the unique vacuum state having $`p=0`$.

## 4 Computer Automation

When $`s=2`$, i.e. in 1-dimensional space and 1-dimensional time, a system of light-cone coordinates is appropriate . In this 2-dimensional space-time, computer automation has been accomplished as deterministic exact analog computation  (computation over “cells” in the continuum of $`^n`$) to obtain primitive extended tube domains of holomorphy for $`W(n;z)`$. By a series of abstractions the computation is done with essentially reversible logic, programming in the Prolog language, and simulating on a Turing machine.

Just as the classical computer, Turing machine, computes over $``$ or equivalently over $`_2`$, we now have what can be called a complex Turing machine, in fact a, severally complex Turing machine.

The primitive extended tube domains are bounded by analytic hypersurfaces, namely Riemann cuts denoted by $`C_{ij}`$ and other hypersurfaces denoted by $`S_{ij,kl}`$ and $`F_{ij,kl}`$ . These domains are in the form of semi-algebraic sets. Since the computation is symbolic, it is also exact, which is important in handling analytic functions.

Because of Lorentz invariance properties of the physics involved, the domains have a structure referred to as complex Lorentz projective spaces 777This is different from Euclidean complex projective spaces, well known in mathematics. The difference is captured in the Hall-Wightman theorem .. Related to this invariance are certain continuum cells over which the computation occurs. Thus this computation is also like analog computation which would otherwise be regarded as impossible to do exactly.

## 5 Analytic Extensions

In relativistic quantum field theory it is possible to implement the physical requirement of microcausality. There exists quantum microcausality (field operators commute or anticommute) at totally space-like points.

Together with the requirement of permutation invariance of the domains, the edge-of-the-wedge theorem provides enlargements of the original primitive domains of analyticity to analyticity into unions of permuted primitive domains.

Mapping these union domains creates some Boolean satisfiability problems. In fact, the novel methods of computation raise interesting issues of computability and complexity.

## 6 Non-deterministic Holomorphic Extensions

By the nature of analytic domains in more than one complex variable, it is in general possible to further extend these domains towards the maximal domains called envelopes of holomorphy. By considering boundary related semi-algebraic sets, there are non-deterministic computations of holomorphic extensions of domains. After the guessing step, the verification is by deterministic processes mentioned above.

Built-in permutation invariance has considerable power, just as $`n!`$ rapidly dominates over $`2^n`$ for large $`n`$.

## 7 Uniformity of Computation

Uniformity in the direction of universal computation has been discussed in , in different contexts, including numerical analysis. We do have certain types of uniformity here.

First we note that the computation is independent of any particular form of Lagrangian or dynamics, and is uniform in $`n`$, qualifying for a universal quantum machine over $`^{\mathrm{}}`$. The latter space is the infinite discrete union $`_{n=1}^{\mathrm{}}^n`$.

### 7.1 Function Order Uniformity

When the progrm runs for $`s=2`$, dynamic memory allocation is used through the operating system. Because $`n`$ can be input as a variable, only part of the whole memory management cost is outside the program. The program itself is independent of $`n=sm`$ and therefore is uniform in $`n`$, which is unbounded above. We can call this function order uniformity in $`n_{\mathrm{}}`$.

### 7.2 Space-time Uniformity

In addition, there is uniformity in the dimension $`s2`$ of space-time, in the following manner. Given a dimension $`s2`$ of space-time, looking at the semi-algebraic sets defining the primitive extended tube domains of holomorphy (hypersurface boundaries), and at function orders, there are three different classes of orders. These classes comprise, a) lower order W functions, b) intermediate order W functions, and c) high order W functions . Extended tube domains for all high order W functions have the same complicacy. For a) we have $`ms+1`$, and for c), $`m>s(s1)/2+2`$. The remaining cases lie in class b). For example, there is no class b) for $`s=2`$, the most complicated primitive domain being for the 3-point function. If $`s=3`$, then $`m=5`$ is the only case in class b). When $`s=4`$, we have in class b), the cases, $`m=6,7`$ and $`8`$.

Since $`s2`$ is unbounded above, we can call this space-time dimension uniformity in $`s_{\mathrm{}}`$.

### 7.3 co-NHolo Uniformity

The holomorphy envelopes $`H[D_m]`$ for different orders $`m`$ of Wightman functions are related  in the following way.

For $`0<r<m`$, and relative to $`H[D_m]`$,

$$H[D_m]\underset{\sigma \mathrm{P}_m}{}\left\{H[\sigma D_{mr}]\times (\sigma ^{sr})\right\},$$

where $`\sigma `$ denotes permutations in $`\mathrm{P}_m`$, the permutation group in the $`m`$ points of the $`m`$-point W function. In the case of Schlicht domains (analogous to single sheeted Riemann surfaces in $``$), the $``$ sign means set theoretic inclusion.

For example, in $`s=2`$, the 4-point function cannot be continued beyond the 2-point function Riemann cuts nor the (permuted) 3-point function Källén-Wightman domains of holomorphy.

This is a statement regarding analyticity that does not exist, and thus refers to the complements of the domains of holomorphy; hence the use of the prefix co-. Because computations of analytic extensions of domains are non-deterministic (hence the notation N), we can say that we have co-NHolo uniformity.

## 8 Conclusion

We started with relativity and continuum quantum field theory. Arbitrary numbers of particles (optionally with spins) can be created or annihilated. Relying on a fruitful set of models, we have related what appeared to be different models of quantum and classical computation based on non-relativistic quantum mechanics and classical mechanics. Exact deterministic and non-deterministic computation over continuous domains appear naturally. Furthermore there is uniformity in computation over, unbounded above, or arbitrarily high order $`n`$ of $`W(n;z)`$ and arbitrarily high dimension $`s`$ of space-time. In the present context this can be called co-NHolo uniformity in $`n_{\mathrm{}}`$ and $`s_{\mathrm{}}`$. The novel methods of computation raise interesting issues of computability and complexity, and possibly could shed more light on quantum field theory itself.


text:
# 1 Chiral order parameter: 𝐵(𝑝̃_𝑘) obtained from Eq. (), as a function of (𝑇,𝜇): (𝑇_𝑐=158 MeV,𝜇=0), (𝑇=0,𝜇_𝑐=275 MeV). As described in the text, the phase boundary is fixed by the condition ℬ(𝑇,𝜇)=0.

TEMPERATURE, CHEMICAL POTENTIAL
AND THE $`\rho `$–MESON
C.D. ROBERTS and S.M. SCHMIDT

Physics Division, Argonne National Laboratory,
Argonne IL 60439-4843, USA

1. Introduction. Models of QCD must confront nonperturbative phenomena such as confinement, dynamical chiral symmetry breaking (DCSB) and the formation of bound states. In addition, a unified approach should describe the deconfinement and chiral symmetry restoring phase transition exhibited by strongly-interacting matter under extreme conditions of temperature and density. Nonperturbative Dyson-Schwinger equation (DSE) models  provide insight into a wide range of zero temperature hadronic phenomena; e.g., non-hadronic electroweak interactions of light- and heavy-mesons , and diverse meson-meson  and meson-nucleon  form factors. This is the foundation for their application at nonzero-($`T,\mu `$,-. Herein we describe the calculation of the deconfinement and chiral symmetry restoring phase boundary, and the medium dependence of $`\rho `$-meson properties. We also introduce an extension to describe the time-evolution in the plasma of the quark’s scalar and vector self energies based on a Vlasov equation.

2. Dyson-Schwinger Equation at Nonzero-($`T,\mu `$). The dressed-quark DSE at nonzero-($`T,\mu `$) is

$`S^1(\stackrel{~}{p}_k)=i\stackrel{}{\gamma }\stackrel{}{p}A(\stackrel{~}{p}_k)+i\gamma _4\omega _{k_+}C(\stackrel{~}{p}_k)+B(\stackrel{~}{p}_k)=i\stackrel{}{\gamma }\stackrel{}{p}+i\gamma _4\omega _{k_+}+\mathrm{\Sigma }(\stackrel{~}{p}_k),`$ (1)
where $`\stackrel{~}{p}_k=(\stackrel{}{p},\omega _{k_+})`$, $`\omega _{k_+}=\omega _k+i\mu `$, and $`\omega _k=(2k+1)\pi T`$ is the quark’s Matsubara frequency. The complex-valued scalar functions: $`A(\stackrel{}{p},\omega _{k_+})`$, $`B(\stackrel{}{p},\omega _{k_+})`$ and $`C(\stackrel{}{p},\omega _{k_+})`$, depend only on $`(|\stackrel{}{p}|^2,\omega _{k_+}^2)`$. With a given dressed-gluon propagator the solutions are determined by

$`B(\stackrel{~}{p}_k)m_0`$ $`=`$ $`{\displaystyle \frac{8}{3}}{\displaystyle \frac{d^4\stackrel{}{q}}{(2\pi )^4}D(\stackrel{~}{p}_k\stackrel{~}{q}_k)\frac{B(\stackrel{~}{q}_k)}{\stackrel{~}{q}_k^2C^2(\stackrel{~}{q}_k)+B^2(\stackrel{~}{q}_k)}},`$ (2)
$`(C(\stackrel{~}{p}_k)1)\stackrel{~}{p}_k^2`$ $`=`$ $`{\displaystyle \frac{4}{3}}{\displaystyle \frac{d^4\stackrel{}{q}}{(2\pi )^4}D(\stackrel{~}{p}_k\stackrel{~}{q}_k)\frac{\stackrel{~}{p}_k\stackrel{~}{q}_kC(\stackrel{~}{q}_k)}{\stackrel{~}{q}_k^2C^2(\stackrel{~}{q}_k)+B^2(\stackrel{~}{q}_k)}},`$ (3)
where herein we only consider models where $`A(p)=C(p)`$. It is the interplay between the functions $`B`$ and $`C`$ that leads to confinement, realised via the absence of a Lehmann representation for the dressed-quark $`2`$-point function ,-. $`B0`$ in the chiral limit signals DCSB.

To provide an illustrative solution of the quark DSE we employ an Ansatz for the scalar function characterising the dressed-gluon propagator :

$$D(p)=3\pi ^2\frac{\eta ^2}{T}\delta _{k0}\delta ^3(p).$$
(4)
The infrared enhancement in this choice ensures quark confinement and DCSB. As an infrared dominant model, Eq. (4) does not represent well the interaction away from $`(\stackrel{~}{p}_k\stackrel{~}{q}_k)^20`$ and that introduces some model-dependent artefacts. However, they are easily identified and the model yields qualitatively reliable results, preserving features of more sophisticated studies.

Using Eq. (4) in Eqs. (2-3) we obtain a system with two phases. The Nambu-Goldstone (NG) phase is characterised by dynamically broken chiral symmetry and confined dressed-quarks. The alternative Wigner-Weyl (WW) solution describes a phase of the model with restored chiral symmetry and deconfinement. In studying the phase transition one must consider the relative stability of the confined and deconfined phases, which is measured by the $`(T,\mu )`$-dependent pressure difference between the two distinct phases: $`(T,\mu )=P[S_{\mathrm{NG}}]P[S_{\mathrm{WW}}]`$. $`(T,\mu )>0`$ indicates the stability of the confined (Nambu-Goldstone) phase and hence the phase boundary is specified by that curve in the $`(T,\mu )`$-plane for which $`(T,\mu )=0`$. The critical line is depicted in Fig. 1. The phase transition is first order for any non-zero $`\mu `$ and second order for $`\mu =0`$. The model has mean field critical exponents, which is a feature of the rainbow-ladder truncation . The study of thermodynamic properties shows that it is essential to keep scalar and vector self-energies as well as their momentum dependence .

Mesons are quark-antiquark bound states and their masses are obtained by solving the Bethe-Salpeter equation . Here we focus on the vector channel and employing Eq. (4) the eigenvalue equation for the bound state mass is

$$\frac{\eta ^2}{2}\mathrm{𝖱𝖾}\left\{\sigma _B(\omega _{0+}^2\frac{1}{4}M_{\rho \pm }^2)^2\left[\pm \omega _{0+}^2\frac{1}{4}M_{\rho \pm }^2\right]\sigma _C(\omega _{0+}^2\frac{1}{4}M_{\rho \pm }^2)^2\right\}=1,$$
(5)
where $`\sigma _{B,C}(\stackrel{~}{p}_k^2)=\{B(\stackrel{~}{p}_k^2),C(\stackrel{~}{p}_k^2)\}/[\stackrel{~}{p}_k^2C^2(\stackrel{~}{p}_k^2)+B^2(\stackrel{~}{p}_k^2)].`$ The equation for the $`\rho `$-meson’s transverse component is obtained with $`[\omega _{0+}^2{\displaystyle \frac{1}{4}}M_\rho ^2]`$ in Eq. (5) and in the chiral-limit yields $`M_\rho ^2={\displaystyle \frac{1}{2}}\eta ^2,\text{independent}\text{ of }T\text{ and }\mu \text{.}`$ This is the $`T=0=\mu `$ result of Ref. . Even for nonzero current-quark mass, $`M_\rho `$ changes by less than 1% as $`T`$ and $`\mu `$ are increased from zero toward their critical values. Its insensitivity is consistent with the absence of a constant mass-shift in the transverse polarization tensor for a gauge-boson. For the longitudinal component one obtains in the chiral limit:

$$M_{\rho +}^2=\frac{1}{2}\eta ^24(\mu ^2\pi ^2T^2).$$
(6)
The results for the medium-dependence of the $`\rho `$ meson are summarised in Fig. 2. As in the case of the dressed-quark mass function, the response to increasing $`T`$ and $`\mu `$ is anti-correlated: the $`\rho `$\- mass decreases with increasing chemical potential and increases with temperature. This anti-correlation leads to an edge along which the $`T`$ and $`\mu `$ effects compensate and the mass remains unchanged up to the transition point.

3. Nonequilibrium Application. The time evolution of the self energies can be studied using Vlasov’s equation

$$_tf(p,x)+_pE(p,x)_xf(p,x)_xE(p,x)_pf(p,x)=0.$$
(7)
Solving this equation is complicated for two reasons. (i) The energy is a functional of the scalar and vector self energies, which in general are nonzero and momentum-dependent. While the scalar self energy is small in the plasma phase due to chiral symmetry restoration, the vector self energy remains significant . (ii) The absence of a Lehmann representation for the dressed-quark propagator in the confined phase precludes the existence of a single particle distribution function, $`f`$, in this phase. Therefore a conventional kinetic theory is only reasonable in the deconfined phase. This situation is adequately represented in DSE models; e.g., Refs.  describe a quark’s $`(T,\mu )`$-evolution from a confined to a propagating mode, and Ref.  makes use of this evolved quasiparticle behaviour in calculating the plasma’s thermodynamic properties. Therefore, approaching the phase boundary from the plasma domain we anticipate a discontinuous disappearance of the quark distribution function, $`f`$.

As an illustration we employ an instantaneous interaction of the form

$$D(p)=3\pi ^2\eta \delta ^3(p)$$
(8)
to represent dynamics in the deconfined phase. In this case the Matsubara sum in Eq. (2) can be performed analytically and we obtain:

$`\mathrm{\Sigma }^B(p,x)`$ $`=`$ $`\eta {\displaystyle \frac{\mathrm{\Sigma }^B(p,x)+m_0}{(1+\mathrm{\Sigma }^C(p,x))E^{}(p,x)}}[12f(p,x)],`$ (9)
$`\mathrm{\Sigma }^C(p,x)`$ $`=`$ $`\eta {\displaystyle \frac{1}{(1+\mathrm{\Sigma }^C(p,x))E^{}(p,x)}}[12f(p,x)],`$ (10)
with the quasi particle energy: $`E^{}(p,x)=\sqrt{(\stackrel{}{p}^{})^2+M^{}(p,x)^2}`$, the renormalized momentum: $`\stackrel{}{p}^{}=\stackrel{}{p}(1+\mathrm{\Sigma }^C(p,x))`$, and mass: $`M^{}(p,x)=m_0+\mathrm{\Sigma }^B(p,x)`$. As a test whether this simplification still yields necessary and qualitatively important features, such as $`C1,Bm_0`$, in Fig. 3 we compare the momentum dependence obtained in the models specified by Eqs. (4,8) in the vicinity of $`T_c`$. Both functions are well reproduced and hence Eq. (8) can be used to model the persistence of non-perturbative effects in the deconfined domain. The solution of Eqs. (7,9-10) provide the time-evolution of the quark self-energy and distribution function.

As in the case of thermal equilibrium, the vector self energy plays an important role. Neglecting $`\mathrm{\Sigma }^C`$ and the momentum dependence of $`\mathrm{\Sigma }^B`$ a simpler equation is obtained

$$_tf(p,x)+\frac{\stackrel{}{p}}{E(p,x)}_xf(p,x)m(x)_xm(x)_pf(p,x)=0,$$
(11)
with $`m(x)`$ the quark mass obtained as a solution of the gap equation in models without confinement. This equation has been widely studied; e.g. Refs. . However, we anticipate that the numerical solution of Eq. (7) will yield significantly different results because of the presence and persistence of the vector self energy in the deconfined domain.

Acknowledgments. This work was supported by the US Department of Energy, Nuclear Physics Division, under contract no. W-31-109-ENG-38, the National Science Foundation under grant nos. INT-9603385 and PHY97-22429, and benefited from the resources of the National Energy Research Scientific Computing Center. S.M.S. acknowledges financial support from the A. v. Humboldt foundation.


text:
# Energy level dynamics in systems with weakly multifractal eigenstates: equivalence to 1D correlated fermions at low temperatures

## Abstract

It is shown that the parametric spectral statistics in the critical random matrix ensemble with multifractal eigenvector statistics are identical to the statistics of correlated 1D fermions at finite temperatures. For weak multifractality the effective temperature of fictitious 1D fermions is proportional to $`T_{eff}(1d_n)/n1`$, where $`d_n`$ is the fractal dimension found from the $`n`$-th moment of inverse participation ratio. For large energy and parameter separations the fictitious fermions are described by the Luttinger liquid model which follows from the Calogero-Sutherland model. The low-temperature asymptotic form of the two-point equal-parameter spectral correlation function is found for all energy separations and its relevance for the low temperature equal-time density correlations in the Calogero-Sutherland model is conjectured.

The spectral statistics in complex quantum systems are signatures of the underlying dynamics of the corresponding classical counterpart. The spectral statistics in chaotic and disordered systems in the limit of infinite dimensionless conductance $`g`$ is described by the classical random matrix theory of Wigner and Dyson  (WD statistics). The WD statistics possess a remarkable property of universality: it depends only on the symmetry class with respect to the time-reversal transformation $`𝒯`$. The three symmetry classes correspond to the lack of $`𝒯`$-invariance (the unitary ensemble, $`\beta =2`$); the $`𝒯`$-invariant systems with $`𝒯^2=1`$ (the orthogonal ensemble, $`\beta =1`$), and the $`𝒯`$-invariant systems with $`𝒯^2=1`$ (the symplectic ensemble, $`\beta =4`$), respectively. The physical ground behind this universality is the structureless eigenfunctions in the ergodic regime which implies the invariance of the eigenfunction statistics with respect to a unitary transformation of the basis.

In real disordered metals the eigenfunctions are not basis-invariant. The basis-preference reaches its extreme form for the strongly impure metals where all eigenfunctions are localized in the coordinate space but delocalized in the momentum space. In this case the spectral statistics is Poissonian in the thermodynamic (TD) limit.

For low-dimensional systems $`d=1,2`$, where all states are localized in the TD limit, one can observe the smooth crossover from the WD to the Poisson spectral statistics as a function of the parameter $`\xi /L`$, where $`\xi `$ is the localization radius and $`L`$ is the system size. The dependence of the spectral correlation functions on the energy variable $`s=E/\mathrm{\Delta }`$ ($`\mathrm{\Delta }`$ is the mean level separation) is non-universal for finite $`L/\xi `$ but all of them tend to the Poisson limit as $`L/\xi \mathrm{}`$.

In systems of higher dimensionality $`d>2`$ the situation is different because of the presence of the Anderson localization transition at a critical disorder $`W=W_c`$. In the metal phase $`W<W_c`$ the dimensional conductance $`g(L)\mathrm{}`$ as $`L\mathrm{}`$ and one obtains the WD spectral statistics in the TD limit. In the insulator state $`g(L)0`$ at $`L\mathrm{}`$ and the limiting statistics is Poissonian. However, there is a fixed point $`W=W_c`$ in which the spectral statistics are nearly independent of $`L`$. Thus at the critical point there exist universal spectral statistics which are neither WD nor Poissonian but rather a hybrid of both . However, the universality of the critical spectral statistics (CSS) is somewhat limited, since it depends not only on the Dyson symmetry parameter $`\beta `$ but also on the critical value of the dimensional conductance $`g^{}`$ which in turn depends on the dimensionality $`d`$ of the system . Thus for each universality class there is a family of critical spectral statistics parametrized by the critical dimensionless conductance $`g^{}`$.

The very existence of the subject of critical level statistics imposes a constraint on the possible values of the localization length exponent $`\nu =[d\beta (g)/d\mathrm{ln}g]^1|_{g=g^{}}`$, where $`\beta (g)d\mathrm{ln}g/d\mathrm{ln}L`$ is the scaling function. Indeed, for the spectral statistics to be meaningful the width of the critical energy window $`\delta E`$ must be much larger than the mean level separation $`\mathrm{\Delta }1/L^d`$. The quantity $`\delta E`$ is defined as the distance from the mobility edge $`E=E_c`$ at which the localization or correlation radius $`\xi (\delta E)|\delta E|^\nu `$ is equal to the system size $`L`$. The number of critical eigenstates $`𝒩=\delta E/\mathrm{\Delta }`$ is proportional to $`L^{d\frac{1}{\nu }}`$. For $`\nu d>1`$ this number tends to infinity in the TD limit $`L\mathrm{}`$ despite the width of the critical energy window shrinks to zero. This necessary condition for the existence of the critical statistics is secured by the famous Harris criterion $`\nu d>2`$.

However, the critical exponent $`\nu `$ enters not only in the necessary condition for the CSS but also in the correlation functions of the density of energy levels $`\rho (E)`$. It has been shown in  that there is a power-law tail in the critical two-level correlation function (TLCF) $`R(\omega )=\rho (E)\rho (E+\omega )`$ that arises because of the finite-size correction to the dimensionless conductance $`g(L_\omega )/g^{}1(L_\omega /L)^{1/\nu }=s^{\frac{1}{\nu d}}`$, where $`L_\omega \omega ^{1/d}L`$ is the length scale set by the energy difference $`\omega =EE^{}s\mathrm{\Delta }`$ between two levels. The sign of this tail depends on whether the critical energy levels are on the metal ($`E<E_c`$) or on the insulator ($`E>E_c`$) side of the mobility edge, in the same way as for the 2D systems in the weak delocalization ($`\beta =4`$) or weak localization ($`\beta =1`$) regimes . Clearly, this power-law tail does not reflect properties of the critical eigenstates but rather the behavior $`\xi (|EE_c|)`$ of the size of the space region where eigenstates show the critical space correlations.

In order to study the relationship between the properties of the critical eigenstates and the CSS in its pure form one should consider a system with a continuous line of critical points where $`\beta (g)=0`$. This case formally corresponds to $`\nu =\mathrm{}`$ and the finite-size effects are absent.

Another complication which makes ambiguous the definition of CSS is the fact that being independent of the system size, the spectral correlation functions depend on the boundary conditions  and topology of a system. Therefore we will consider the system of the torus topology where CSS takes its ‘canonical’ form. In particular, the TLCF decays exponentially in this case.

As has been already mentioned the universality of the WD statistics is based on the ergodic, basis-invariant statistics of eigenfunctions which one may encounter in different physical situations. The characteristic feature of all critical quantum systems is the multifractal statistics of the critical eigenfunctions . The simplest two-point correlations of the critical wave functions can be obtained from the renormalization group result  for $`l<r<\xi <L`$:

$$|\mathrm{\Psi }_E(0)|^{2n}|\mathrm{\Psi }_E(r)|^{2n}=p<|\mathrm{\Psi }_E(0)|^{2n}>^2(\xi /r)^{\alpha _n},$$
(1)
where $`l`$ is the short-distance cut-off of the order of the elastic scattering length, $`p=(\xi /L)^d`$ is the probability for a reference point to be inside a localization region. The exponent $`\alpha _n=2n(d_nd_{2n})+d_{2n}2d_n+d`$ is expressed through the fractal dimensions $`d_n`$ defined by the $`L`$-dependence of the moments of inverse participation ratio:

$$L^d|\mathrm{\Psi }_E(r)|^{2n}L^{d_n(n1)},n2.$$
(2)
At the critical point $`p=1`$ and the correlation radius $`\xi \mathrm{}`$ in Eq.(1) must be replaced by the sample size $`L`$. Spectral statistics are related with eigenfunction correlations at different energies $`E`$ and $`E^{}=E+\omega `$. If $`\omega \mathrm{\Delta }`$ and thus $`L_\omega L`$ one should substitute $`L_\omega `$ for $`L`$ in the $`r`$-dependent term of Eq.(1). In this way the multifractality exponents $`d_n`$ enter the spectral $`\omega `$-dependences .

For weak multifractality one can expect the fractal dimensions to be a linear function of $`n`$ which is controlled by only one parameter $`a`$:

$$d_n/d=1an,a1/g^{}1.$$
(3)
This relationship holds approximately for the Anderson transition in $`2+ϵ`$ dimensions  and is fulfilled exactly for the critical eigenstate of the Dirac equation in the random vector-potential .

Thus for critical quantum systems with weak multifractality it is natural to expect that the spectral statistics depends on only one system-specific parameter - the critical conductance $`g^{}`$.

In view of the expected universality, it is useful to find a simple one-parameter random matrix ensemble with the multifractal eigenfunction statistics which would play the same role for the critical systems as the classical RMT does for the ergodic systems. As a matter of fact there are few candidates . However, here we focus only on one of them , since for this ensemble the multifractality of eigenstates has been rigorously proven .

Consider a Hermitean $`N\times N`$ matrix with the real ($`\beta =1`$), complex ($`\beta =2`$) or quaternionic ($`\beta =4`$) entries $`H_{ij}`$ $`(ij)`$ which are independent Gaussian random numbers with zero mean and the variance:

$$|H_{ij}|^2=\frac{1}{1+\frac{(ij)^2}{B^2}}\{\begin{array}{cc}1/\beta ,& i=j,\\ 1/2,& ij\end{array}J(ij)$$
(4)
This model has been shown to be critical both for large  and for small  values of $`B`$ with the fractal dimensionality $`d_2`$ at the center of the spectral band being:

$$d_2=\{\begin{array}{cc}1\frac{1}{\pi \beta B},& B1,\\ 2B,& B1\end{array}$$
(5)
Thus the 1D system with long-range hopping described by the matrix Hamiltonian Eq.(4) possesses the line of critical points $`B(0,\mathrm{})`$, the fractal dimensionality $`d_2`$ changing from 1 to 0 with decreasing $`B`$.

One can extend this matrix 1D model by closing it into a ring and applying a ‘flux’ $`\phi [0,1]`$. In this case

$$H_{ij}(\phi )=H_{ij}+H_{ij}^{(1)}e^{2\pi i\phi sgn(ij)}$$
(6)
is a sum of two independent Gaussian random numbers with the variance of $`H_{ij}^{(1)}`$ given by:

$$|H_{ij}^{(1)}|^2=J(N|ij|).$$
(7)
For large values of $`B`$ which correspond to weak multifractality one can derive  an effective field theory – the supersymmetric nonlinear sigma-model  – which describes the spectral and eigenfunction correlations of the critical random matrix ensemble Eq.(4):

$$F[𝐐]=\frac{g^{}}{16}\underset{i,j=1}{\overset{N}{}}Str[𝐐_iU_{|ij|}𝐐_j]+\frac{i\pi s}{4N}\underset{i=1}{\overset{N}{}}Str[\sigma _z𝐐_i],$$
(8)
where $`𝐐`$ is the supermatrix with $`𝐐_i^2=\mathrm{𝟏}`$ and

$$g^{}=4\beta B.$$
(9)
The symmetry with respect to time reversal is encoded in the symmetry of $`𝐐_i`$ in exactly the same way as for the diffusive sigma-model . The only difference is the long-range kernel $`U_{|ij|}`$ with the Fourier-transform $`\stackrel{~}{U}_k=|k|`$. For a torus geometry $`k=2\pi m/N`$, where $`m`$ is an arbitrary integer.

One can explicitly resolve the constraint $`𝐐^2=1`$ by switching to the integration over the ‘angles’ $`W`$. Then the Gaussian fluctuations of ‘angles’ recover the spectrum of ‘quasi-diffusion’ modes:

$$\epsilon _m=g^{}|m|,m=0,\pm 1,\pm 2,\mathrm{}$$
(10)
The problem of spectral statistics can be generalized to include the dependence of spectrum on the flux $`\phi `$ introduced by Eq.(6). One can define  the parametric two-level correlation function $`R(s,\phi )=\rho (E,0)\rho (E+s\mathrm{\Delta },\phi )`$ which can be treated in the framework of the same nonlinear sigma-model but with the phase-dependent $`\epsilon _m`$:

$$\epsilon _m(\phi )=g^{}|m\phi |.$$
(11)
Following the work by Andreev and Altshuler Ref. we introduce the spectral determinant:

$$D^1(s,\phi )=\underset{m0}{}\frac{\epsilon _m^2(\phi )+s^2}{\epsilon _m^2(0)}$$
(12)
Then it can be shown in the same way as in Ref. that the parametric TLCF for $`s1`$ and $`g^{}\phi 1`$ can be expressed in terms of the spectral determinant as follows:

$`R^u(s,\phi )`$ $`=`$ $`{\displaystyle \frac{1}{4\pi ^2}}{\displaystyle \frac{^2G(s,\phi )}{s^2}}+\mathrm{cos}(2\pi s)e^{G(s,\phi )}`$ (13)
$`R^o(s,\phi )`$ $`=`$ $`{\displaystyle \frac{1}{2\pi ^2}}{\displaystyle \frac{^2G(s,\phi )}{s^2}}+2\mathrm{cos}(2\pi s)e^{2G(s,\phi )}`$ (14)
$`R^s(s,\phi )`$ $`=`$ $`{\displaystyle \frac{1}{8\pi ^2}}{\displaystyle \frac{^2G(s,\phi )}{s^2}}+{\displaystyle \frac{\pi }{\sqrt{8}}}\mathrm{cos}(2\pi s)e^{G(s,\phi )/2}+`$ (15)
$`+`$ $`{\displaystyle \frac{1}{8}}\mathrm{cos}(4\pi s)e^{2G(s,\phi )},`$ (16)
where $`G(s,\phi )=G(s,\phi +1)`$ is a periodic in $`\phi `$ function:

$$e^{G(s,\phi )}=\frac{D(s,\phi )}{2\pi ^2(s^2+\epsilon _0^2(\phi ))}.$$
(17)
Eqs.(13-16) coincide with the corresponding formulae in Ref. for the unitary, orthogonal and symplectic ensembles after some misprints are corrected as in Ref. and $`s2s`$ for the symplectic ensemble to take account of the Kramers degeneracy. The only difference is in the form of the spectral determinant Eq.(12) due to the specific spectrum $`\epsilon _m(\phi )`$ of the quasi-diffusion modes.

Using the functional representation Eq.(8) one can also find the leading term in the deviation from the WD statistics at $`sg^{}`$ and $`\phi =0`$ using the results of Ref.:

$$\delta R(s)=\frac{1}{2\pi ^2\beta }\left(\underset{m0}{}\frac{1}{\epsilon _m^2}\right)\frac{d^2}{ds^2}\left[s^2R_{WD}(s)\right],$$
(18)
where $`R_{WD}(s)`$ is the Wigner-Dyson TLCF.

A remarkable property of the function $`G(s,\phi )`$ for the critical sigma-model Eq.(8) on a torus is that it can be decomposed into the sum $`G(s,\phi )=F(z)+F(\overline{z})`$ of analytic functions $`F(z)`$ and $`F(\overline{z})`$ where $`\tau =g^{}\phi `$, $`z=\tau +is`$, $`\overline{z}=\tau is`$ and

$$F(z)=\mathrm{ln}[\sqrt{2}g^{}\mathrm{sin}(\pi z/g^{})].$$
(19)
Eq.(19) results from a straightforward evaluation  of the product in Eq.(12).

On the other hand, it can be easily verified that $`G(s,\phi )`$ given by Eq.(19) is proportional to the Green’s function of the free-boson field $`\mathrm{\Phi }(z)`$ on a torus in the $`(1+1)`$ $`z`$-space: $`0<\mathrm{}z<g^{}`$, $`\mathrm{}<\mathrm{}z<+\mathrm{}`$:

$$\mathrm{\Phi }(s,g^{}\phi )\mathrm{\Phi }(0,0)_S\mathrm{\Phi }(0,0)\mathrm{\Phi }(0,0)_S=KG(s,\phi ),$$
(20)
where $`\mathrm{}_S`$ denotes the functional average with the free-boson action:

$$S[\mathrm{\Phi }]=\frac{1}{8\pi K}_0^g^{}𝑑\tau _{\mathrm{}}^+\mathrm{}𝑑s\left[(_s\mathrm{\Phi })^2+(_\tau \mathrm{\Phi })^2\right].$$
(21)
Now we are in a position to make a crucial step and suggest that for the critical RMT described by Eq.(4), the Andreev-Altshuler equations (13-16) are nothing but density-density correlations in the Luttinger liquid of fictitious 1D fermions at a finite temperature $`T=1/g^{}`$:

$$R(s,\phi )=\overline{n}^2n(s,\tau )n(0,0)_S1,\tau =g^{}\phi .$$
(22)
Indeed, the density operator $`n(s,\tau )`$ ($`s`$\- is space and $`\tau (0,1/T)`$ is imaginary time coordinate) for 1D interacting fermions with the Fermi-momentum $`k_F=\pi `$ can be expressed through the free boson field $`\mathrm{\Phi }(s,\tau )`$ as follows :

$`n(s,\tau )`$ $`=`$ $`{\displaystyle \frac{1}{2\pi }}_s\mathrm{\Phi }(s,\tau )+A_K\mathrm{cos}[2\pi s+\mathrm{\Phi }(s,\tau )]+`$ (23)
$`+`$ $`B_K\mathrm{cos}[4\pi s+2\mathrm{\Phi }(s,\tau )].`$ (24)
The constants $`A_K`$ and $`B_K`$ are independent of ‘temperature’ $`1/g^{}`$ but depend on the interaction constant $`K`$. They can be uniquely determined from the WD limit $`g^{}\mathrm{}`$.

Using Eqs.(22,23) and the well known result for the Gaussian average of the exponent:

$$e^{ip\mathrm{\Phi }(s,\tau )}e^{ip\mathrm{\Phi }(0,0)}=e^{Kp^2[G(s,\tau )G(0,0)]},$$
(25)
one can verify that for the choice:

$$K=\frac{2}{\beta },\beta =1,2,4$$
(26)
the Andreev-Altshuler formulae Eqs.(13-16) are reproduced exactly for the orthogonal, unitary and symplectic ensembles, respectively. Now we remind on a known result that the parametric spectral statistics in the WD limit $`g^{}=\mathrm{}`$ is equivalent to the Tomonaga-Luttinger liquid at zero temperature. It follows directly from Ref. and the equivalence of the Calogero-Sutherland model and the Tomonaga-Luttinger liquid for large distances $`|z|1`$. The critical random matrix ensemble Eq.(4,6,7) and the critical 1D sigma-model Eq.(8) turns out to be the simplest generalization of the WD theory that retains the Tomonaga- Luttinger liquid analogy extended for finite ‘temperatures’ $`T=1/g^{}`$ which are related with the spectrum of fractal dimensions Eq.(3).

The Wigner-Dyson two-level statistics for all three symmetry classes can be expressed through the single kernel $`K(s)=\mathrm{sin}(\pi s)/(\pi s)`$ in the following way :

$`R^u(s)=K^2(s),`$ (27)
$`R^o(s)=K^2(s){\displaystyle \frac{dK(s)}{ds}}{\displaystyle _s^{\mathrm{}}}K(x)𝑑x,`$ (28)
$`R^s(s)=K^2(s)+{\displaystyle \frac{dK(s)}{ds}}{\displaystyle _0^s}K(x)𝑑x.`$ (29)
It turns out that such a representation is also valid for the critical TLCF at $`g^{}1`$ if the kernel is replaced by:

$$K(s)=\frac{T}{\alpha }\frac{\mathrm{sin}(\pi \alpha s)}{\mathrm{sinh}(\pi Ts)},T=\frac{1}{g^{}}.$$
(30)
where $`\alpha =1`$ for the orthogonal and the unitary ensemble and $`\alpha =2`$ for the symplectic ensemble. The form of the kernel Eq.(28) can be guessed from the well known density correlation function for the case of free fermions in one dimension at a finite temperature that corresponds to the unitary ensemble . In order to prove this statement we note that for $`s1`$ and $`g^{}1`$ Eqs.(27-29) with the kernel Eq.(30) give the same leading terms as Eqs.(13-16) with $`G(s,0)`$ given by Eq.(19). For $`sg^{}`$ Eqs.(27-29) give corrections to the WD statistics that coincide with Eq.(18). Thus the representation Eqs.(27-29) with the kernel Eq.(30) are the correct asymptotic expressions for both $`s1`$ and $`sg^{}`$. At $`g^{}1`$ these regions have a parametrically large overlap so that Eqs.(27-29) are valid for all values of $`s`$. Given that for all $`s`$ at $`T=0`$ and for large energy separations $`s1`$ at $`T1`$ expressions Eqs.(27-29) correspond to equal-time correlations in the Calogero-Sutherland model, one can expect Eqs.(27-29) with the kernel Eq.(30) to describe the low-temperature equal-time correlations in the Calogero-Sutherland model  for for all $`s`$ and the interaction constant $`K=2/\beta `$.

One can use the kernel Eq.(28) to compute the spacing distribution function $`P(s,g^{})`$ which can be compared with the corresponding distribution function $`P_c(s)`$ obtained by the numerical diagonalization of the three-dimensional Anderson model at the critical point. Such a comparison is done in Ref. for $`\beta =1,2,4`$. It turns out that identifying the parameter $`g^{}`$ from the fitting $`P(s,g^{})e^{\kappa (g^{})s}`$ with the far exponential tail of $`P_c(s)e^{\kappa s}`$ one reproduces the entire distribution function $`P_c(s)`$ extremely well.

A. M. T. acknowledges a kind hospitality of Abdus Salam ICTP where this work was performed.


text:
# Finite-genus solutions for the Hirota’s bilinear difference equation.

## Abstract

The finite-genus solutions for the Hirota’s bilinear difference equation are constructed using the Fay’s identities for the $`\theta `$-functions of compact Riemann surfaces.

In the present work I want to consider once more the question of constructing the finite-genus solutions for the famous Hirota’s bilinear difference equation (HBDE)

$$\tau (k1,l,m)\tau (k+1,l,m)+\tau (k,l1,m)\tau (k,l+1,m)+\tau (k,l,m1)\tau (k,l,m+1)=0$$
(1)
which has been solved in  using the so-called algebraic-geometrical approach. This method, which is the most powerful method for deriving the quasi-periodic solutions (QPS) and which has been developed for almost all known integrable systems, exploits some rather sophisticated pieces of the theory of functions of complex variables and is based on some theorems determining the number of functions with prescribed structure of singularities on the Riemann surfaces (see  for review). However, in some cases the QPS can be found with less efforts, in a more straightforward way, using the fact that the finite-genus QPS (and namely they are the subject of this note) of all integrable equations possess similar and rather simple structure: up to some phases they are ’meromorphic’ combinations of the $`\theta `$-functions associated with compact Riemann surfaces  (the situation resembles the pure soliton case, where solutions are rational functions of exponents). Thus, to construct these solutions we only have to determine some constant parameters. This, as in the pure soliton case, can be done directly, using the well-known properties of the $`\theta `$-functions of compact Riemann surfaces. In  such approach was developed for the Ablowitz-Ladik hierarchy, where the finite-genus solutions were ’extracted’ from the so-called Fay’s formulae . The fact, that the HBDE is closely related to the Fay’s identities is not new and was mentioned, e.g., in  (see Remark 2.7), but contrary to this work I will use these identities as a starting point and will show how one can derive from them, by rather short and very simple calculations (which to my knowledge have not been presented explicitly in the literature), a wide family of solutions for the HBDE.

Consider a compact Riemann surface $`X`$ of the genus $`g`$. One can choose a set of closed contours (cycles) $`\{a_i,b_i\}_{i=1,\mathrm{},g}`$ with the intersection indices

$$a_ia_j=b_ib_j=0,a_ib_j=\delta _{ij}i,j=1,\mathrm{},g$$
(2)
and find $`g`$ independent holomorphic differentials which satisfy the normalization conditions

$$_{a_i}\omega _k=\delta _{ik}$$
(3)
The matrix of the $`b`$-periods,

$$\mathrm{\Omega }_{ik}=_{b_i}\omega _k$$
(4)
determines the so-called period lattice, $`L_\mathrm{\Omega }=\left\{𝒎+\mathrm{\Omega }𝒏,𝒎,𝒏^g\right\}`$, the Jacobian of this surface $`\mathrm{Jac}(X)=^g/L_\mathrm{\Omega }`$ (2$`g`$ torus) and the Abel mapping $`X\mathrm{Jac}(X)`$,

$$𝓐(P)=_{P_0}^P𝝎$$
(5)
where $`𝝎`$ is the $`g`$-vector of the 1-forms, $`𝝎=(\omega _1,\mathrm{},\omega _g)^T`$ and $`P_0`$ is some fixed point of $`X`$.

A central object of the theory of the compact Riemann surfaces is the $`\theta `$-function, $`\theta (𝜻)=\theta (𝜻,\mathrm{\Omega })`$,

$$\theta \left(𝜻\right)=\underset{𝒏^g}{}\mathrm{exp}\left\{\pi i(𝒏,\mathrm{\Omega }𝒏)+2\pi i(𝒏,𝜻)\right\}$$
(6)
where $`(𝝃,𝜼)`$ stands for $`_{i=1}^g\xi _i\eta _i`$, which is a quasiperiodic function on $`^g`$

$`\theta \left(𝜻+𝒏\right)`$ $`=`$ $`\theta \left(𝜻\right)`$ (7)
$`\theta \left(𝜻+\mathrm{\Omega }𝒏\right)`$ $`=`$ $`\mathrm{exp}\left\{\pi i(𝒏,\mathrm{\Omega }𝒏)2\pi i(𝒏,𝜻)\right\}\theta \left(𝜻\right)`$ (8)
for any $`𝒏^g`$.

The famous Fay’s trisecant formula can be written as

$$\underset{i=1}{\overset{3}{}}a_i\theta \left(𝜻+𝜼_i\right)\theta \left(𝜻𝜼_i\right)=0$$
(9)
where

$`2𝜼_1`$ $`=`$ $`𝓐(P_1)+𝓐(P_2)+𝓐(P_3)𝓐(P_4),`$ (10)
$`2𝜼_2`$ $`=`$ $`𝓐(P_1)𝓐(P_2)+𝓐(P_3)𝓐(P_4),`$ (11)
$`2𝜼_3`$ $`=`$ $`𝓐(P_1)+𝓐(P_2)𝓐(P_3)𝓐(P_4).`$ (12)
Here $`P_1,\mathrm{},P_4`$ are arbitrary points of $`X`$, and the constants $`a_i`$ are given by

$`a_1`$ $`=`$ $`e(P_4,P_1)e(P_2,P_3),`$ (13)
$`a_2`$ $`=`$ $`e(P_4,P_2)e(P_3,P_1),`$ (14)
$`a_1`$ $`=`$ $`e(P_4,P_3)e(P_1,P_2).`$ (15)
The skew-symmetric function $`e(P,Q)`$, $`e(P,Q)=e(Q,P)`$, is closely related to the prime form  and is given by

$$e(P,Q)=\theta [𝜹^{},𝜹^{\prime \prime }]\left(𝓐(Q)𝓐(P)\right)$$
(16)
where $`\theta [𝜶,𝜷]\left(𝜻\right)`$ is the so-called $`\theta `$-function with characteristics,

$$\theta [𝜶,𝜷]\left(𝜻\right)=\mathrm{exp}\left\{\pi i(𝜶,\mathrm{\Omega }𝜶)+2\pi i(𝜶,𝜻+𝜷)\right\}\theta \left(𝜻+\mathrm{\Omega }𝜶+𝜷\right),$$
(17)
and $`(𝜹^{},𝜹^{\prime \prime })\frac{1}{2}^{2g}/^{2g}`$ is a non-singular odd characteristics,

$$\theta [𝜹^{},𝜹^{\prime \prime }]\left(\mathrm{𝟎}\right)=0,\mathrm{grad}_𝜻\theta [𝜹^{},𝜹^{\prime \prime }]\left(\mathrm{𝟎}\right)\mathrm{𝟎}$$
(18)
Now it is very easy to establish relations between (9) and the HBDE (1). To do this one has first to introduce the discrete variables by

$$\mathrm{\Theta }(k,l,m)=\theta \left(𝜻+k𝜼_1+l𝜼_2+m𝜼_3\right)$$
(19)
The Fay’s identity (9) can now be rewritten as

$`a_1\mathrm{\Theta }(k1,l,m)\mathrm{\Theta }(k+1,l,m)+a_2\mathrm{\Theta }(k,l1,m)\mathrm{\Theta }(k,l+1,m)`$ (20)
$`+a_3\mathrm{\Theta }(k,l,m1)\mathrm{\Theta }(k,l,m+1)`$ $`=`$ $`0`$ (21)
from which it follows that the quantity

$$\tau (k,l,m)=a_1^{k^2/2}a_2^{l^2/2}a_3^{m^2/2}\mathrm{\Theta }(k,l,m)$$
(22)
satisfies HBDE (1).

Thus the last formula,

$$\tau (k,l,m)=a_1^{k^2/2}a_2^{l^2/2}a_3^{m^2/2}\theta \left(𝜻+k𝜼_1+l𝜼_2+m𝜼_3\right)$$
(23)
where $`\theta `$ is the $`\theta `$-function of some compact Riemann surface $`X`$ and the constants $`a_i`$ and $`\eta _i`$ depend on four points (paths) on this surface, determines a family of finite-genus solutions of the HBDE.

Now I want to discuss the solutions obtained above. First of all it should be noted that these solutions are ’finite-genus’ but not quasiperiodic. At first glance, this seems to be strange, because usually the finite-genus solutions naturally appear when one solves quasiperiodic problems. For, example, in $`1+1`$ dimensional discrete systems, such as Toda chain, Ablowitz-Ladik equations, etc, the quasiperiodicity leads to the polynomial dependence of the scattering matrix of the auxiliary problem on the spectral parameter. These polynomials determine some hyperelliptic curve (spectral curve) of finite genus, and the quasiperiodic solutions are built up of the $`\theta `$-functions corresponding to the latter. Thus, in some sense, in discrete systems the quasiperiodicity implies the ’finite-genus’ property. In other words, quasiperiodic solutions are finite-genus. However, the reverse is not obligatory true. The Fay’s identities do not imply that the integrals $`𝓐(P)`$ are in some way related to the periods $`𝝎`$. Formula (9) determine, so to say, ’local’ properties of the $`\theta `$-functions, and all above consideration was local, without reference to some boundary conditions (quasiperiodicity).

Another point which I would like to discuss here is to compare the approach of this note with the algebro-geometrical one. The ideology of the latter is to operate on the Riemann surface: the key moment in calculating the Baker-Akhiezer function, the central object of the algebro-geometrical method, is to satisfy the condition that it is a single-valued function of the point of the Riemann surface. Since we didn’t introduce the Baker-Akhiezer function and didn’t study its analytical properties the question of how our solutions depend on the points $`P_iX`$ (or on the integral paths from $`P_0`$ to $`P_i`$, to be more precise) is not so crucial as in the algebro-geometrical approach and is in some sense the question of parametrization of constants. One can restrict the points $`P_i`$, $`i=1,2,3`$, to some neighborhood of the point $`P_0`$ and rewrite the Abel’s integrals in terms of local coordinates (with some polynomial representing the Riemann surface). This enables not to mention the Riemann surface and reformulate all results in terms of integrals over the complex plane. As to the ’global’ (or ’homotopical’) effects, which arise when we add to the paths $`(P_0,P_i)`$ some integer cycles ($`_{k=1}^gm_ka_k+n_kb_k`$, $`m_k,n_k`$) it should be noted that, if we consider the Abel integral as mapping $`X^g`$, then such deformations of the contours change $`𝓐(P)`$ as $`𝓐(P)𝓐(P)+𝜸`$, $`𝜸L_\mathrm{\Omega }`$. This results first in adding some half-period to the argument of the $`\theta `$-function in (23),

$$𝜻+\underset{i=1}{\overset{3}{}}k_i𝜼_i𝜻+\underset{i=1}{\overset{3}{}}k_i𝜼_i+\frac{1}{2}\underset{i=1}{\overset{3}{}}k_i𝜸_𝒊,𝜸_𝒊L_\mathrm{\Omega }$$
(24)
and, second, in altering the constants $`a_i`$’s, $`a_ia_i\mathrm{exp}(f_i)`$. Thus, we come to the point where one can try to apply the theory of Backlund transformations for the HBDE to describe the transformations of $`\theta `$-functions due to half-period shifts (and vice versa). This is an interesting problem, which deserves special studies.

To conclude I would like to note the following. In principle, the direct approach based on the Fay’s identities can be used to derive the finite-genus solutions not only for the HBDE but for almost all known integrable systems (some examples one can find in the book ). However, contrary to the case of the HBDE where all ’calculations’ take only two lines (formulae (19) and (22) above), in the case of partial differential equations such as, e.g., KP equation the corresponding calculations become rather cumbersome. For example, to solve the KP equation one has to expand the Fay’s identities up to the third order in some small parameter. From the other hand, it is a widely known fact that almost all known integrable equations can be derived from the HBDE. Hence one can ’skip’ the Fay’s identities and use solutions (23) as a starting point. In  it was shown how to obtain some finite-genus solutions for the nonlinear Schrödinger and KP equations using the corresponding solutions for the Ablowitz-Ladik hierarchy, which can be viewed as some ’pre-continuous’ version of the HBDE (in the Ablowitz-Ladik hierarchy 2 of 3 discrete coordinates of the HBDE are presented as the Miwa’s shifts of two infinite sets of continuous variables).


text:
# Induced instability for boson-fermion mixed condensate of Alkali atoms due to attractive boson-fermion interaction

## A Life Time of Collective Tunneling Effect in Metastable Region

To estimate a life time of metastable mixed condensates, we consider the collective tunneling effect with the collective variable $`R`$ (the boson radius) and its effective energy $`E(R)`$ as a collective potential. The effective energy $`E(R)`$ is obtained from (14), and, for the metastable condensates, a typical shape of it is the line (c) in Fig. 1. We approximate this $`E(R)`$ with the linear-plus-harmonic-oscillator type potential eq. (16) (Fig. 3). $`M`$ is inertia mass for the collective variable $`R`$. The meanings of other parameters ($`F`$, $`G`$, $`\mathrm{\Omega }`$, $`R_{eq}`$, $`R_t`$) can be read off in Fig. 3, and they are fixed in order that the potential $`V(R)`$ reproduce the $`E(R)`$ numerically obtained by (14). This approximation should be enough for the order estimation of the life time. For the kinetic term, we take

$$T=\frac{\mathrm{}^2}{2M}\frac{d^2}{dR^2}.$$
(A1)
In $`V(R)`$, the metastable state $`\psi _M(x)`$ before tunneling is approximately given by the ground-state wave function in the harmonic oscillator potential $`V_2(R)`$:

$$\psi _M(R)=\frac{1}{\sqrt{\pi ^{\frac{1}{2}}a_{HO}}}\mathrm{exp}\left[\frac{(RR_{eq})^2}{a_{HO}^2}\right],$$
(A2)
where $`a_{HO}=\sqrt{\mathrm{}/M\mathrm{\Omega }}`$ is a harmonic oscillator length. This state has the zero-point energy $`\mathrm{}\mathrm{\Omega }/2`$ measured from $`V_0`$. The continuum state $`\psi _D(R)`$ after the tunneling decay is obtained by the wave function in the linear potential $`V_1(R)`$, and its Schrödinger equation becomes

$$\left[\frac{\mathrm{}^2}{2M}\frac{d^2}{dR^2}+E(FRG)\right]\psi _D(R)=0.$$
(A3)
To solve the above equation, we use the WKB approximation. It should be noticed that, because the state energy is lower than the potential maximum ($`E<FR_tG`$), a turning point exist at $`R_E=(E+G)/F`$. Thus, the WKB connection formula should be used at $`R_E`$ for the continuum state $`\psi _D`$:

$`\psi _D(R)`$ $`=`$ $`{\displaystyle \frac{A}{2\sqrt{K(R)}}}\mathrm{sin}\left[{\displaystyle \frac{2}{3}}\left({\displaystyle \frac{2MF}{\mathrm{}^2}}\right)^{\frac{1}{2}}(R_ER)^{\frac{3}{2}}+{\displaystyle \frac{\pi }{4}}\right],\text{for }0R<R_E`$ (A4)
$`\psi _D(R)`$ $`=`$ $`{\displaystyle \frac{A}{2\sqrt{K(R)}}}\mathrm{exp}\left[{\displaystyle \frac{2}{3}}\left({\displaystyle \frac{2MF}{\mathrm{}^2}}\right)^{\frac{1}{2}}(RR_E)^{\frac{3}{2}}\right],R_E<R`$ (A5)
where $`A=\sqrt{2M/(\pi \mathrm{}^2)}`$ and $`\mathrm{}K(R)=\sqrt{2MF|R_ER|}`$.

Let us consider the tunneling decay rate $`\mathrm{\Gamma }`$ from $`\psi _M(R)`$ to $`\psi _D(R)`$. Using the golden rule of Fermi, it becomes

$$\mathrm{\Gamma }=\frac{2\pi }{\mathrm{}}\left|𝑑R\psi _D^{}(R)[V(R)V_2(R)]\psi _M(R)\right|_{E=E_M}^2,$$
(A6)
where $`E_M=V_0+\mathrm{}\mathrm{\Omega }/2`$ is the energy of $`\psi _M(R)`$. Because $`V(R)V_2(R)=0`$ for $`R_t<R`$ and $`\psi _M(0),\psi _M^{}(0)1`$, eq. (A6) becomes

$$\mathrm{\Gamma }=\frac{2\pi }{\mathrm{}}\left(\frac{\mathrm{}^2}{2M}\right)^2\left[\frac{d\psi _D(R)}{dR}\psi _M(R)\psi _D(R)\frac{d\psi _M(R)}{dR}\right]_{R=R_t}.$$
(A7)
Using eqs. (A2, A4, A5), the tunneling rate $`\mathrm{\Gamma }`$ is shown to have the form $`\mathrm{\Gamma }=D\mathrm{exp}(W)`$ : $`D`$ can be interpreted as the staying probability Per unit time and $`\mathrm{exp}(W)`$ is the transition coefficient. The explicit formulae for $`D`$ and $`W`$ are

$`D`$ $`=`$ $`{\displaystyle \frac{4\sqrt{2}\mathrm{\Omega }}{\sqrt{\pi }}}\sqrt{{\displaystyle \frac{\mathrm{\Delta }V}{\mathrm{}\mathrm{\Omega }}}{\displaystyle \frac{1}{2}}}`$ (A8)
$`W`$ $`=`$ $`{\displaystyle \frac{4}{3}}\sqrt{{\displaystyle \frac{2M}{\mathrm{}^2}}}\left(\mathrm{\Delta }V{\displaystyle \frac{\mathrm{}\mathrm{\Omega }}{2}}\right)(R_tR_E)+{\displaystyle \frac{(R_0R_t)^2}{a_{HO}^2}},`$ (A9)
where $`\mathrm{\Delta }V=V_1(R_t)V_0`$ is a barrier height. To represent the parameters $`F`$ and $`G`$ in $`V_1`$ by $`\mathrm{\Delta }V`$, we have used

$$V_1(R_t)=FR_tG,V_1(R_E)=FR_EG=E=V_0+\frac{1}{2}\mathrm{}\mathrm{\Omega }.$$
(A10)
The tunneling life-time is defined as an inverse of decay rate: $`\tau _{ct}=\mathrm{\Gamma }^1`$.

## B Critical Condition for Unstable Region

In this appendix, we rederive the Mølmer scaling law which gives the critical condition for the unstable region . We assume the Thomas-Fermi approximations both for the boson and fermion density distributions. In that approximations, the density distributions are given as solutions of the Thomas-Fermi equations:

$`\stackrel{~}{g}n_b(x)+x^2+\stackrel{~}{h}n_f(x)=\stackrel{~}{\mu }_b,`$ (B1)
$`[6\pi ^2n_f(x)]^{\frac{2}{3}}+x^2+\stackrel{~}{h}n_b(x)=\stackrel{~}{\mu }_f,`$ (B2)
where $`n_b(x)`$ is the boson density distribution scaled by the harmonic-oscillator length $`\xi =(\mathrm{}/m\omega )^{1/2}`$ ($`n_b(x)=|\mathrm{\Phi }(x)|^2\xi ^3`$). The scaled chemical potentials, $`\stackrel{~}{\mu }_{b,f}`$, and coupling constants, $`\stackrel{~}{g}`$ and $`\stackrel{~}{h}`$, have been defined in the main body of this paper.

Eliminating $`n_b(x)`$ in (B2) by (B1), we obtain the equation $`F[n_f(x)]=G[n_f(x)]`$ to determine the fermion density $`n_f(x)`$, where

$$F[n_f(x)][6\pi ^2n_f(x)]^{\frac{2}{3}},G[n_f(x)]\stackrel{~}{\mu }_f+\frac{|\stackrel{~}{h}|}{\stackrel{~}{g}}\stackrel{~}{\mu }_b\left(1+\frac{|\stackrel{~}{h}|}{\stackrel{~}{g}}\right)x^2+\frac{|\stackrel{~}{h}|^2}{\stackrel{~}{g}}n_f(x).$$
(B3)
We concentrate on the central density $`n_f(0)`$. In order to determine $`n_f(0)`$ two conditions should be satisfied:

$$F[n_f(0)]=G[n_f(0)],\frac{\delta F}{\delta n_f}[n_f(0)]\frac{\delta G}{\delta n_f}[n_f(0)],$$
(B4)
Evaluating equations in (B4) with (B3), we obtain the critical condition for the unstable region:

$$\frac{|\stackrel{~}{h}|^2}{\stackrel{~}{g}}=\frac{4\pi ^2}{\sqrt{3}}\left(\stackrel{~}{\mu }_f+\frac{|\stackrel{~}{h}|}{\stackrel{~}{g}}\stackrel{~}{\mu }_b\right).$$
(B5)
To evaluate the right-hand side of (B5), we should use the relations between $`\stackrel{~}{\mu }_{b,f}`$ and the boson/fermion particle number $`N_{b,f}`$, which are obtained by solving the Thomas-Fermi equations in (B2). Here we assume that $`\stackrel{~}{\mu }_b=0`$ and $`\stackrel{~}{\mu }_f=2(6N_f)^{1/3}`$ ($`\stackrel{~}{\mu }_f`$ for a free fermion system) for (B5). Consequently, we obtain the Mølmer scaling relation:

$$\frac{|\stackrel{~}{h}|^2}{\stackrel{~}{g}}=\frac{4\pi ^2}{6^{1/6}\sqrt{6}}N_f^{1/2}12.0N_f^{1/6}.$$
(B6)
It should be noted that the coefficient 12.0 in (B6) are close to the value 13.8 which is obtained by Mølmer  with numerically evaluating the Thomas-Fermi equations.


text:
# TeV Emission by Ultra-High Energy Cosmic Rays in Nearby, Dormant AGNs

The curvature radiation produced by particles accelerating near the event horizon of a spinning supermassive black hole, threaded by externally supported magnetic field lines is considered. It is shown that light nuclei suffer catastrophic curvature losses that limit the maximum energy they can attain to values well below that imposed by the maximum potential difference induced by the black hole dynamo, unless the curvature radius of magnetic field lines largely exceeds the gravitational radius. It is further shown that the dominant fraction of the rotational energy extracted from the black hole is radiated in the TeV band. Given the observed flux of ultra-high energy cosmic rays, and the estimated number of nearby supermassive black holes, it is expected that if dormant AGNs are the sources of UHECRs, as proposed recently by Boldt & Ghosh, then they should also be detectable at TeV energies by present TeV experiments.

Energy losses resulting from interactions with the cosmic microwave background limit the distance a cosmic ray (CR) of energy $`>10^{20}`$ eV can traverse to less than 50 Mpc . As a consequence, if the origin of the ultra-high energy cosmic rays (UHECRs) observed is associated with astrophysical objects, rather than decaying supermassive X particles as in the top-down scenario , then the sources of UHECRs must be close by. Plausible classes of UHECR sources have been discussed in the literature, including active galactic nuclei (AGNs)  and gamma-ray bursts .

Recently, it has been proposed  that the UHECR observed at Earth may originate from dormant AGNs in the local Universe having masses in excess of $`10^9M_{}`$. The idea being that these objects, although underluminous relative to active quasars, contain spinning supermassive black holes that, instead of producing luminous radio jets as seen in blazars, serve as accelerators of a small number of particles to ultra-high energies. The acceleration mechanism invoked is associated with a BZ type process ; specifically, individual particles are accelerated by the potential difference induced by the rotation of a black hole that is threaded by externally produced magnetic field lines, during episodes when breakdown of the vacuum does not occur. In the absence of energy losses, the maximum energy a particle can gain by this mechanism is limited to the voltage drop involved, and is proportional to the product of the black hole mass and the strength of magnetic field.

A recent analysis  indicates that massive dark objects (MDOs) are present in the centers of nearby galaxies, some of which have masses in excess of $`10^{10}M_{}`$. A plausible interpretation is that the MDOs are supermassive black holes, and may represent quasar remnants or dormant AGNs. This interpretation is further supported by the fact that a correlation between the black hole mass and bulge luminosity, similar to that found for the sample of galaxies studied in ref. , has been found  for a sample of bright quasars, using a completely different method. By applying their model to a list of objects from ref. , of which 14 have compact central masses larger than $`10^9M_{}`$, Boldt & Ghosh  estimate that protons can be accelerated up to energies in excess of a few times $`10^{20}`$ eV. As emphasized by them, in order to account for the measured flux of UHECRs, an average power of only $`10^{42}`$ erg s-1 is needed to be extracted from a black hole by the accelerated cosmic rays, corresponding to a mass loss rate of order $`10^{10}`$ g s-1, and a loss rate of electric charge that constitutes only a small fraction of the total effective current required to induce the potential difference across the gap.

Acceleration to the maximum energy allowed is possible provided that radiative losses are sufficiently small. Boldt & Ghosh  argued that proton energy losses during the acceleration phase due to pair production and photomeson production on ambient photons are unimportant by virtue of the low accretion luminosity anticipated in those objects. However, they have not discussed the radiation associated with the acceleration process itself. In the following, we consider the curvature radiation produced by the accelerated particles, and show that in the case of light nuclei it limits the maximum energy attainable to values below the full voltage, unless the average radius of curvature of the particle’s trajectory exceeds the gravitational radius by at least a factor of a few. We further show that the curvature radiation is emitted predominantly in the TeV band, and conclude that if dormant AGNs are indeed the sources of UHECRs, then they should also be detectable by present, ground based experiments at TeV energies.

Curvature radiation: The electric potential difference generated by a maximally rotating black hole of mass $`M=10^9M_9M_{}`$, threaded by magnetic field of strength $`B=10^4B_4`$ Gauss is

$$\mathrm{\Delta }V4.4\times 10^{20}B_4M_9(h/R_g)^2\mathrm{volts},$$
(1)
where $`h`$ is the gap height, and $`R_g=GM/c^2`$ is the gravitational radius. In the presence of a nonuniform magnetic field, particles accelerated by this potential difference will suffer energy losses through curvature radiation, even if initially they move along magnetic field lines. Since the gyroradius of a proton having energy $`ϵ`$, $`R_c=ϵ/eB`$, is smaller than the gravitational radius:

$$R_c/R_g=(ϵ/e\mathrm{\Delta }V)e\mathrm{\Delta }Vc^2/(eGBM)(ϵ/e\mathrm{\Delta }V)<1,$$
(2)
(the gyroradius of an ion of charge $`Z`$ having energy near the maximum imposed by the voltage drop will be larger by a factor of $`Z`$), we expect the average radius of curvature of a particle’s trajectory to be of order the curvature radius of magnetic field lines in the gap. (The curvature radii of different trajectories should span some range though, reflecting the different boundary conditions.) The computations of particles’ trajectories, even for relatively simple magnetic field topologies, are complicated by the fact that the gyroradius at the highest energies is comparable to the size of the hole, as can be seen from eq. (2), and are beyond the scope of this paper. In what follows, we denote by $`\rho `$ the average curvature radius of an accelerating ion, and assume that it is independent of the ion energy. The rate of energy loss through curvature radiation by a particle of energy $`ϵ=mc^2\gamma `$ can then be expressed as

$$P=\frac{2}{3}\frac{e^2c\gamma ^4}{\rho ^2}.$$
(3)
The energy change per unit length of an accelerating ion having charge $`Z`$ and mass $`m_i=\mu m_p`$ is given by

$$dϵ/ds=eZ\mathrm{\Delta }V/hP/c,$$
(4)
yielding a maximum acceleration energy,

$$ϵ_{max}=3\times 10^{19}\mu Z^{1/4}M_9^{1/2}B_4^{1/4}(\rho ^2h/R_g^3)^{1/4}\mathrm{eV},$$
(5)
where eq. (3) has been used. Consequently, only a fraction

$$\eta =0.1\mu M_9^{1/2}(ZB_4)^{3/4}(\rho /R_g)^{1/2}(h/R_g)^{7/4};\eta 1,$$
(6)
of the potential energy available will be released as UHECRs; the rest will be radiated in the form of curvature photons. For the most massive MDOs listed in table 2 of ref.  ($`M_9>20`$) we obtain $`ϵ_{max}1.5\times 10^{20}\mu (ZB_4\rho ^2h/R_g^3)^{1/4}`$ eV, and $`\eta 0.02\mu (ZB_4)^{3/4}(\rho /R_g)^{1/2}(h/R_g)^{7/4}`$. Thus $`B_4\rho ^2h/R_g^320`$ is required in order to accelerate a proton ($`Z=\mu =1`$) to energies $`3\times 10^{20}`$ eV in these systems, corresponding to $`\eta =10^3(\rho /R_g)^2(h/R_g)^1`$. The requirement on heavier nuclei is more relaxed. Estimates of the maximum value of the horizon threading magnetic field strength yield $`B_41`$ , and according to recent numerical simulations $`B_4`$ may be well below unity . Consequently, acceleration of light nuclei to the highest energies measured by current experiments requires $`\rho `$ to be larger than $`R_g`$ by a factor of at least a few (assuming $`hR_g`$). The values of $`B_4`$ obtained by Boldt & Ghosh  for the Magorian et al. sample , assuming equipartition between the magnetic field in the vicinity of the horizon and the matter infalling into the center, lie in the range between 0.1 to 1, given the estimated mass loss rate of the galaxy bulge.

The spectrum produced by the curvature radiation of a single ion will peak at an energy

$$ϵ_{\gamma max}=1.5\gamma ^3\mathrm{}c/\rho =1.6\times 10^7ϵ_{max}\mu ^1(ZB_4)^{1/2}(h/R_g)^{1/2}=5M_9^{1/2}(ZB_4)^{3/4}(\rho ^2h^3/R_g^5)^{1/4}\mathrm{TeV},$$
(7)
and is a power law $`I(ϵ_\gamma )(ϵ_\gamma )^{1/3}`$ below the cutoff. The overall spectrum of curvature photons would depend on the energy distribution of the accelerating particles, and is expected to be somewhat softer below the peak. For $`ϵ_{max}=3\times 10^{20}`$ eV and $`hR_g`$ we obtain from eq. (7), $`ϵ_{\gamma max}50\mu ^1(ZB_4)^{1/2}`$ TeV.

The number of TeV photons per proton produced in the process is roughly

$$n_\gamma eZ\mathrm{\Delta }V/ϵ_{\gamma max}=10^8M_9^{1/2}(ZB_4)^{1/4}(h/R_g)^{5/4}(\rho /R_g)^{1/2}.$$
(8)
The mean free path to pair creation of a photon moving at an angle $`\chi `$ to the magnetic field is $`l95(B_4\mathrm{sin}\chi )^1e^q`$, with $`q=1.3\times 10^3M_9^{1/2}B_4^{7/4}\mathrm{sin}\chi ^1(\rho ^2h^3/R_g^5)^{1/4}`$ and, therefore, the radiated photons escape the system freely for the range of parameters considered here. However, when $`q`$ becomes sufficiently small (at $`q`$ of about 50 one photon per proton will be converted into an electron-positron pair), a pair cascade may be initiated with high enough probability to lead to a breakdown of the gap.

Observational consequences: As shown recently , the observed CR spectrum above $`10^{19}`$ eV, as measured by the AGASA  and Fly’s Eye  experiments, can be accounted for by a homogeneous cosmological distribution of CR sources with power law spectra and energy production rate of $`1.5\times 10^{37}`$ ergs Mpc-3 s-1. Given this energy production rate, and denoting by $`n_{CR}`$ the density of objects contributing to the observed CR flux in this energy range, the average power released in the form of UHECR by a single source can be expressed as

$$L_{CR}=2\times 10^{42}\left(\frac{n_{CR}}{10^4Mpc^3}\right)^1\mathrm{erg}\mathrm{s}^1.$$
(9)
Employing eq. (9), one finds that the average TeV flux emitted by a single CR source at a distance of $`D=50D_{50}`$ Mpc is

$$F_\gamma 10^{12}\left(\frac{n_{CR}}{10^4Mpc^3}\right)^1\eta ^1D_{50}^2\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1,$$
(10)
where $`\eta <1`$ is the UHECR production efficiency defined in eq. (6). As seen from eq. (6), $`\eta `$ should be of order unity if the particles accelerated by the black hole dynamo are predominantly heavy nuclei. There is, however, evidence that the CR composition changes from heavy nuclei below the ankle (at energy of $`5\times 10^{18}`$ eV) to light nuclei above it . Assuming a protonic CR composition above the ankle, we obtain values of $`\eta `$ between 10-1 and 10-3 for the range of parameters considered above.

The density of black holes above a certain mass can be estimated using the correlation between bulge luminosity and MDO mass found by Magorian et al. , and an appropriate luminosity function of nearby galaxies to correct for the incompleteness of their sample. However, as seen from eq. (5), the maximum energy a particle can attain depends on a combination of parameters and not solely on the mass. Therefore, the black hole mass above which particles can accelerate to the required energies is uncertain, and since the MDOs having relevant masses lie on the bright end of the luminosity function, this uncertainty in mass translates to a large uncertainty in $`n_{CR}`$. For a reasonable choice of parameters we anticipate that $`M_9`$ of at least a few would be required in order to account for UHECR energies observed. Using the luminosity function measured by Efstathiou et al. , we estimate that the density of MDOs having masses $`M_9>1`$ is of order a few times $`10^4`$ Mpc-3. An estimate based on a k-band luminosity function measured more recently  yielded a similar value.

The threshold flux for a 5$`\sigma `$ detection of gamma-rays above TeV by current TeV experiments is $`5\times 10^{12}t_{day}^{1/2}`$ erg s-1 cm -2, where $`t_{day}`$ is the exposure time measured in days . With the above estimation of $`\eta `$ and $`n_{CR}`$, we expect that at least some fraction of the UHECR sources will be detectable by present TeV experiments, provided that the TeV photons escape the system. Conceivable sources of opacity that may give rise to attenuation of the TeV flux are considered next.

The curvature photons produced near the black hole can be absorbed through pair production on IR photons in the galaxy. The corresponding optical depth at a given energy depends on the spectrum of IR photons and the energy dependence of the cross section. To an order of magnitude it is given by, $`\tau _{\gamma \gamma }\sigma _{PP}n_{IR}R`$, where $`R`$ is the size of IR emission region - of order several kpc as inferred from the light profiles, $`n_{IR}`$ is the number density of IR photons, and $`\sigma _{PP}`$ is the pair production cross section, given approximately by $`\sigma _{PP}0.2\sigma _T`$ at energies near the threshold. Using this approximation one finds that TeV photons would escape the galaxy provided that the IR luminosity at energy $`ϵ_{IR}=(m_ec^2)^2/ϵ_\gamma `$ satisfies:

$$L_{IR}(ϵ_\gamma ^1)<3\times 10^{45}\left(\frac{R}{3\mathrm{kpc}}\right)\left(\frac{ϵ_\gamma }{1\mathrm{TeV}}\right)^1\mathrm{erg}\mathrm{s}^1.$$
(11)
The integrated IR luminosity implied by eq. (11) would presumably be larger. Condition (11) is satisfied in most cases. A central continuum source, if present, would also contribute a pair production opacity. The IR emission may arise from a cold accretion disk or dust reprocession, and is likely to originate from small radii, of order 10 to 100 $`R_g`$ in these low luminosity objects. Adopting a size of 10 times the gravitational radius of a $`10^{10}M_{}`$ black hole for the IR emission region, we find that the TeV flux will be strongly attenuated if the corresponding IR luminosity exceeds $`10^{40}(ϵ_\gamma /1\mathrm{TeV})^1`$ ergs s-1. This is slightly below the luminosity inferred for low-luminosity AGNs , and comparable to upper limits on the luminosity of the point source in elliptical galaxies .

In conclusion, it has been shown that particles accelerating near the horizon of a spinning suppermassive black hole threaded by externally supported magnetic field lines, suffer severe energy losses through curvature emission. The curvature losses limit the maximum energy attainable by light nuclei to values well below that imposed by the voltage drop. The major fraction of the energy extracted from the rotating hole is radiated in the TeV band, with a rather hard spectrum that extends well beyond 10 TeV. Given the energy flux of cosmic rays above $`10^{19}`$ eV, as reported recently by current CR experiments, and an estimate of the density of supermassive black holes in the universe, it is concluded that if dormant AGNs are the sources of the ultra-high energy cosmic rays, then they should be detectable by current TeV experiments.

I thank Ari Laor, Dan Maoz, Eli Waxman, David Eichler, Amiel Sternberg, and Hagai Netzer for discussions and useful comments. Support from the Israel Science Foundation is acknowledged.