text:
# A particular matrix and its relationship with Fibonacci numbers

## 1 Introduction

Let us consider the $`n\times n`$ symmetric matrix $`𝐀_n`$

$$𝐀_n=\left[\begin{array}{ccccc}1& 1& 1& \mathrm{}& 1\\ 1& 2& 2& \mathrm{}& 2\\ 1& 2& 3& \mathrm{}& 3\\ \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}\\ 1& 2& 3& \mathrm{}& n\end{array}\right],$$

that is

$$a_{ij}=\mathrm{min}(i,j),i,j=1,\mathrm{\hspace{0.17em}2},\mathrm{},n$$

This matrix is, up to a positive scalar, the covariance matrix of a stochastic process with increments which possess the same variance and are uncorrelated.

For future reference consider also the $`(nk+1)\times (nk+1)`$ matrix $`𝐂_{n,k}`$

$$𝐂_{n,k}=\left[\begin{array}{ccccc}k& k& k& \mathrm{}& k\\ k& k+1& k+1& \mathrm{}& k+1\\ k& k+1& k+2& \mathrm{}& k+2\\ \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}\\ k& k+1& k+2& \mathrm{}& n\end{array}\right],k\text{ integer, }1<k<n.$$

For general reference on matrix theory see, for example, , , .

## 2 Determinants

Let us consider the following determinants (dimensions $`n\times n`$)

$$\mathrm{\Delta }_n(i_1,i_2,\mathrm{},i_n)=\left|\begin{array}{ccccc}i_1& i_1& i_1& \mathrm{}& i_1\\ i_1& i_1+i_2& i_1+i_2& \mathrm{}& i_1+i_2\\ i_1& i_1+i_2& i_1+i_2+i_3& \mathrm{}& i_1+i_2+i_3\\ \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}\\ i_1& i_1+i_2& i_1+i_2+i_3& \mathrm{}& i_1+i_2+\mathrm{}+i_n\end{array}\right|,$$

and

$`\mathrm{\Theta }_n(i_1,i_2,\mathrm{},i_n,i_{n+1})=`$
$`=\left|\begin{array}{ccccc}i_1& i_1+i_2& i_1+i_2& \mathrm{}& i_1+i_2\\ i_1& i_1+i_2+i_3& i_1+i_2+i_3& \mathrm{}& i_1+i_2+i_3\\ i_1& i_1+i_2+i_3& i_1+i_2+i_3+i_4& \mathrm{}& i_1+i_2+i_3+i_4\\ \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}\\ i_1& i_1+i_2+i_3& i_1+i_2+i_3+i_4& \mathrm{}& i_1+i_2+\mathrm{}+i_{n+1}\end{array}\right|,`$
where $`i_j,j`$ are integers $`1`$. Let us not that the determinant of $`𝐀_n`$ is $`\mathrm{\Delta }_n(1,\mathrm{\hspace{0.17em}1},\mathrm{},\mathrm{\hspace{0.17em}1})`$ and that of the matrix $`𝐂_{n,k}`$ is $`\mathrm{\Delta }_n(k,\mathrm{\hspace{0.17em}1},\mathrm{},\mathrm{\hspace{0.17em}1})`$.

We are going to prove that

$$\mathrm{\Delta }_n(i_1,i_2,\mathrm{},i_n)=\underset{j=1}{\overset{n}{}}i_j,$$
(2)
and

$$\mathrm{\Theta }_n(i_1,i_2,\mathrm{},i_n,i_{n+1})=i_1\underset{j=3}{\overset{n+1}{}}i_j.$$
(3)
We are going to prove the claim using induction.

First of all

$`\mathrm{\Delta }_2(i_1,i_2)`$ $`=`$ $`\left|\begin{array}{cc}i_1& i_1\\ i_1& i_1+i_2\end{array}\right|`$
$`=`$ $`i_1(i_1+i_2i_1)`$
$`=`$ $`i_1i_2,`$
and

$`\mathrm{\Theta }_2(i_1,i_2,i_3)`$ $`=`$ $`\left|\begin{array}{cc}i_1& i_1+i_2\\ i_1& i_1+i_2+i_3\end{array}\right|`$
$`=`$ $`i_1(i_1+i_2+i_3i_1i_2)`$
$`=`$ $`i_1i_3.`$
Now assume that the claim holds for $`n1`$. If we evaluate $`\mathrm{\Delta }_n(i_1,i_2,\mathrm{},i_n)`$ expanding with respect to the elements of the first row we have that all the algebraic complements except those relative to the first two elements of the row are equal to zero: indeed all the above have identical the first two columns, and we obtain

$`\mathrm{\Delta }_n(i_1,i_2,\mathrm{},i_n)`$ $`=`$ $`i_1\mathrm{\Delta }_{n1}(i_1+i_2i_3,\mathrm{},i_n)i_1\mathrm{\Theta }_{n1}(i_1,i_2,\mathrm{},i_n)`$
$`=`$ $`i_1(i_1+i_2)i_3\mathrm{}i_ni_1i_1i_3\mathrm{}i_n`$
$`=`$ $`i_1i_3i_4\mathrm{}i_n(i_1+i_2i_1)`$
$`=`$ $`{\displaystyle \underset{j=1}{\overset{n}{}}}i_j.`$
Now reverting to $`\mathrm{\Theta }_n(i_1,i_2,\mathrm{},i_n,i_{n+1})`$ the same reasoning as before leads to

$`\mathrm{\Theta }_n(i_1,i_2,\mathrm{},i_n,i_{n+1})`$ $`=`$ $`i_1\mathrm{\Delta }_{n1}(i_1+i_2+i_3i_4,\mathrm{},i_n,i_{n+1})`$
$`(i_1+i_2)\mathrm{\Theta }_{n1}(i_1,i_2+i_3,i_4,\mathrm{},i_{n+1})`$
$`=`$ $`i_1(i_1+i_2+i_3)i_4\mathrm{}i_{n+1}(i_1+i_2)i_1i_4\mathrm{}i_{n+1}`$
$`=`$ $`i_1i_4\mathrm{}i_{n+1}(i_1+i_2+i_3i_1i_2)`$
$`=`$ $`i_1{\displaystyle \underset{j=3}{\overset{n+1}{}}}i_j.`$
This completes the proof.

As a corollary we see that immediately we have $`|𝐀_n|=1`$, and $`|𝐂_{n,k}|=k`$.

## 3 Symmetric Functions of Eigenvalues

Just to fix notation, given $`n`$ numbers

$$\{a_i\}_1^m=\{a_1,a_2,\mathrm{},a_n\}$$

let $`S_k^na_i`$ denote the symmetric functions, that is

$$S_k^na_i=\underset{\stackrel{i_1,\mathrm{},i_k}{i_1<i_2<\mathrm{}<i_k}}{}a_{i_1}a_{i_2}\mathrm{}a_{i_k},1kn.$$

Now given a $`n\times n`$ matrix A with eigenvalues $`\{\lambda _i\}`$ we have

$$S_k^n\lambda _i=\underset{i}{}\left|𝐀_{(i)}^{(nk)}\right|,$$
(6)
where

$$𝐀_{(i)}^{(h)},h=0,\mathrm{\hspace{0.17em}1},\mathrm{},n1$$

denotes one of the principal submatrices of A obtained deleting $`h`$ rows and the corresponding $`h`$ columns. The number of these submatrices is

$$\left(\genfrac{}{}{0pt}{}{n}{h}\right).$$

Note that the submatrices in Equation 6 have dimensions $`k\times k`$.

Considering $`𝐀_n`$ we have

$$S_n^n\lambda _i=|𝐀_n|=1,$$
(7)
$$S_1^n\lambda _i=\mathrm{tr}(𝐀_n)=\frac{n(n+1)}{2}.$$
(8)
A closed inspection of the structure of the determinants in Equation 6 show that

$`S_k^n\lambda _i`$ $`=`$ $`{\displaystyle \underset{i_1=1}{\overset{nk+1}{}}}{\displaystyle \underset{i_2=1}{\overset{nk+2i_1}{}}}{\displaystyle \underset{i_3=1}{\overset{nk+3i_1i_2}{}}}\mathrm{}`$
$`\mathrm{}{\displaystyle \underset{i_k=1}{\overset{ni_1i_2\mathrm{}i_{k1}}{}}}\mathrm{\Delta }_k(i_1,i_2i_3,\mathrm{},i_k)`$
$`=`$ $`{\displaystyle \underset{i_1=1}{\overset{nk+1}{}}}{\displaystyle \underset{i_2=1}{\overset{nk+2i_1}{}}}{\displaystyle \underset{i_3=1}{\overset{nk+3i_1i_2}{}}}\mathrm{}`$
$`\mathrm{}{\displaystyle \underset{i_k=1}{\overset{ni_1i_2\mathrm{}i_{k1}}{}}}i_1i_2\mathrm{}i_k.`$
We can develop this formula according the index $`i_1`$ in the following way

$`S_k^n\lambda _i`$ $`=`$ $`1\times {\displaystyle \underset{i_2=1}{\overset{nk+1}{}}}{\displaystyle \underset{i_3=1}{\overset{nk+2i_2}{}}}\mathrm{}{\displaystyle \underset{i_k=1}{\overset{n1i_2\mathrm{}i_{k1}}{}}}i_2i_3\mathrm{}i_k`$
$`+2\times {\displaystyle \underset{i_2=1}{\overset{nk}{}}}{\displaystyle \underset{i_3=1}{\overset{nk+1i_2}{}}}\mathrm{}{\displaystyle \underset{i_k=1}{\overset{n2i_2\mathrm{}i_{k1}}{}}}i_2i_3\mathrm{}i_k`$
$`+\mathrm{}\mathrm{}`$
$`+(nk+1)\times {\displaystyle \underset{i_2=1}{\overset{1}{}}}{\displaystyle \underset{i_3=1}{\overset{2i_2}{}}}\mathrm{}{\displaystyle \underset{i_k=1}{\overset{k1i_2\mathrm{}i_{k1}}{}}}i_2i_3\mathrm{}i_k.`$
This shows that we can write

$$S_k^n\lambda _i=\underset{i=1}{\overset{nk+1}{}}iS_{k1}^{ni}\lambda _i.$$
(9)
In the above formula we set $`S_0^n\lambda _i=1`$.

Using this equation we have

$$S_k^n\lambda _i=S_{k1}^{n1}\lambda _i+2S_{k1}^{n2}\lambda _i+3S_{k1}^{n3}\lambda _i+\mathrm{}+(nk+1)S_{k1}^{k1}\lambda _i,$$

and

$$S_k^{n1}\lambda _i=S_{k1}^{n2}\lambda _i+2S_{k1}^{n3}\lambda _i+3S_{k1}^{n4}\lambda _i+\mathrm{}+(nk)S_{k1}^{k1}\lambda _i.$$

It follows

$$S_k^n\lambda _i=S_k^{n1}\lambda _i+\underset{i=1}{\overset{nk+1}{}}S_{k1}^{ni}\lambda _i.$$
(10)
This formula makes it easy to evaluate $`S_k^n\lambda _i`$ for all $`n`$ and $`k`$.

Consider the following identity relating to binomial coefficients

$$\left(\genfrac{}{}{0pt}{}{n+k}{nk}\right)=\left(\genfrac{}{}{0pt}{}{n+k1}{nk1}\right)+\underset{i=1}{\overset{nk+1}{}}\left(\genfrac{}{}{0pt}{}{n+k1i}{nk+1i}\right),0kn.$$
(11)
If we set

$$S_k^n=\left(\genfrac{}{}{0pt}{}{n+k}{nk}\right),$$

we see that

$$S_k^{n1}=\left(\genfrac{}{}{0pt}{}{n+k1}{nk1}\right),$$

and

$$S_{k1}^{ni}=\left(\genfrac{}{}{0pt}{}{n+k1i}{nk+1i}\right).$$

In this way Equation 10 follows from identity 11, and conversely. So we can conclude that

$$S_k^n\lambda _i=\left(\genfrac{}{}{0pt}{}{n+k}{nk}\right).$$
(12)
Note that

$$\left(\genfrac{}{}{0pt}{}{n+k}{nk}\right)=\left(\genfrac{}{}{0pt}{}{n+k}{2k}\right).$$

Then Identity 1.76 in  says

$$\underset{k=0}{\overset{n}{}}\left(\genfrac{}{}{0pt}{}{n+k}{nk}\right)=F_{2n+1},$$

where $`F_i`$ is the $`i`$-th Fibonacci number. It follows

$`{\displaystyle \underset{k=1}{\overset{n}{}}}S_k^n\lambda _i+1`$ $`=`$ $`{\displaystyle \underset{k=1}{\overset{n}{}}}\left({\displaystyle \genfrac{}{}{0pt}{}{n+k}{2k}}\right)+1`$ (13)
$`=`$ $`{\displaystyle \underset{k=0}{\overset{n}{}}}\left({\displaystyle \genfrac{}{}{0pt}{}{n+k}{2k}}\right)`$
$`=`$ $`F_{2n+1}.`$
We can get another recurrence relationship. Indeed from

$`{\displaystyle \frac{n+k}{nk}}\left({\displaystyle \genfrac{}{}{0pt}{}{n+k1}{nk1}}\right)`$ $`=`$ $`{\displaystyle \frac{(n+k)(n+k1)!}{(nk)(nk1)!(2k)!}}`$
$`=`$ $`{\displaystyle \frac{(n+k)!}{(2k)!(nk)!}}`$
$`=`$ $`\left({\displaystyle \genfrac{}{}{0pt}{}{n+k}{nk}}\right),`$
we get

$$S_k^n\lambda _i=\frac{n+k}{nk}S_k^{n1}\lambda _i,k=1,\mathrm{\hspace{0.17em}2},\mathrm{},n1$$
(14)


format: text

text:
# Nova V1974 Cygni - results of the 1997 campaign

## 1 Introduction

V1974 Cygni (Nova Cygni 1992) was one of the brightest galactic novae of the 20th century. At the end of February 1992 it reached its maximum brightness equal to $`V=4.3`$ mag.

During subsequent years Nova Cygni 1992 was intensively observed using several techniques. The X-ray flux of the object indicated high luminosity at the level of $`10^{37}`$ erg s-1 caused by hot material burning in thermonuclear reactions at the surface of a white dwarf (Krautter et al. 1996).

The earliest photometric observations obtained by DeYoung and Schmidt (1993, 1994) revealed the presence of a light curve modulation with a period of 0.0813 days (117 min) and an amplitude equal to 0.16 mag.

In 1994, Semeniuk et al. (1994) discovered another periodicity in the light curve of V1974 Cygni with a period equal to 0.085 days (122 min). Subsequent observations showed that this period is not stable and changes with $`\dot{P}`$ at a level of $`10^6`$ (Semeniuk et al. 1995, Olech et al. 1996, Skillman et al. 1997).

The 117 min stable period was interpreted as the orbital period of the system (DeYoung and Schmidt 1994, Semeniuk et al. 1995, Olech et al. 1996). The nature of the 122 min period still remains unexplained. Until now, two major explanations have been considered: 1) it might be either a superhump period of a SU UMa permanent superhumper star (Leibowitz et al. 1995, Skillman et al. 1997, Retter et al. 1997); or 2) it might be either the rotation period of the magnetized white dwarf in the binary system (Semeniuk et al. 1995, Olech et al. 1996).

In this paper we present the results of CCD $`I`$-band photometry of V1974 Cygni performed during the end of the observational season 1996 and for the whole season of 1997.

## 2 Observations

The reported observations of V1974 Cygni were obtained on nine nights during the period from October 4, 1996 to July 15, 1997 at the Ostrowik station of the Warsaw University Observatory.

The 60-cm Cassegrain telescope, equipped with a Tektronics TK512CB back illuminated CCD camera, was used. The scale of the camera was 0.76"/pixel providing a $`6.5^{}\times 6.5^{}`$ field of view. The full description of the telescope and camera was given by Udalski and Pych (1992).

We have monitored the star in the Cousins $`I`$ filter. The exposure times were 180 or 240 seconds depending on the atmospheric conditions, the seeing and the altitude of the object.

A full journal of our CCD observations of V1974 Cygni is given in Table 1.

## 3 Data Reduction

As in our previous papers describing the photometric behavior of V1974 Cygni (Semeniuk et al. 1994, 1995, Olech et al. 1996) the initial reduction of the CCD frames (subtracting BIAS, DARK and FLAT-FIELD images) was done using the standard routines in the IRAF CCDPROC package. 111IRAF is distributed by the National Optical Astronomy Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the National Science Foundation.

A further reduction was done using a newly developed image subtraction method ISIS (Alard and Lupton 1998, Alard 2000). ISIS is based on the fast optimal image subtraction algorithm. Subtraction is done between the reference image (obtained by stacking several best seeing and low background images) and each image from the data set transformed to the coordinate grid of the reference image. Before subtraction the reference image is modified to exactly match the seeing of each image from the data set. This is done by finding the convolution kernel and difference in background levels between the reference and subtracted image. There are more details in the papers of Alard and Lupton (1998) and Alard (2000).

It was shown by Alard (1999) and Olech et al. (1999) that ISIS works very well even in very dense regions of the Galactic bulge and in the cores of the globular clusters. In the first applications ISIS produced light curves with a quality significantly better than obtained by traditional methods such as DaoPHOT and DoPHOT.

During recent years the image subtraction method was widely used for several astronomical objects such as globular clusters (Olech at al. 1999, Kaluzny, Olech and Stanek 2001, Kopacki 2001), open clusters (Mochejska et al. 2002), galaxies (Mochejska et al. 2001), quasars (Woźniak et al. 2000) or microlensing events (Alard 1999, Soszyński et al 2001).

The field of V1974 Cygni located in the Galactic plane is quite rich in stars. There are a few hundreds stellar objects in the $`6.5^{}\times 6.5^{}`$ field of view of the Ostrowik telescope thus CCD frames of the V1974 Cygni are optimal for reduction using the ISIS package.

In the first step, all frames were transformed to a common $`(x,y)`$ coordinate grid. Next we constructed the reference image stacking 14 of the best exposed and low background images. The obtained reference image is shown in the upper panel of Fig. 1.

Image subtraction was then applied to all the frames from our data set. For each CCD frame, the reference image was convolved with a kernel to match its point spread function (PSF), and then it was subtracted from the frame producing set of difference images. One such images is shown in the middle panel of Fig. 1.

At this stage the standard ISIS package produces ABS image, which is built by computing a median of absolute deviations on all subtracted images. Such an image obtained for our data set is shown in the lower panel of Fig. 1. The detection of the potential variables on such a frame can be done using for example the find algorithm of DaoPHOT software (Stetson 1987).

In our case we have already known the position of the variable which we were interested in. In Fig. 1 the variable V1974 Cyg is thus marked by an open circle. One can clearly see that it is detectable by both its difference and ABS image.

The main disadvantage of using the ISIS package is obtaining the light curves in ADU and not in the magnitudes. Thus we should convert the light curves from ADU to magnitudes. For this purpose, the reference frame was reduced with the DaoPHOT/ALLSTAR package (Stetson 1987). Knowing the instrumental magnitude of the variable produced by the ALLSTAR routine and ADU flux of this object at the reference frame, the ADU light curve can be transformed into magnitudes.

The resulting light curve from the night of October 4/5, 1996 is displayed in Fig. 2 and light curves from six of eight July 1997 nights are shown in Fig. 3.

## 4 Power Spectrum Analysis of 1997 Data

To obtain power spectra of the 1997 photometric observations of V1974 Cygni we have used the AoV (Schwarzenberg-Czerny 1989) statistic. Before calculating the power spectrum we have removed the nightly mean and a longer-scale trend from each individual run. Fig. 4 presents results of the analysis for all July 1997 nights.

Upper panel of Fig. 4 shows the AoV power spectrum for the raw light curve. The solid ticks in the figure indicate the frequencies 11.807 and 12.302 cycles/day. The corresponding periods are equal to $`0.08470\pm 0.00014`$ days ($`122.0\pm 0.2`$ min) and $`0.0813\pm 0.0002`$ days ($`117.1\pm 0.3`$ min).

Knowing both periodicities we have fitted our light curve with the Fourier sine series in the following form:

$$mag=A_0+\underset{j=1}{\overset{4}{}}A_j^1\mathrm{sin}(\frac{2\pi jt}{P_1}+\phi _j^1)+\underset{j=1}{\overset{4}{}}A_j^2\mathrm{sin}(\frac{2\pi jt}{P_2}+\phi _j^2)$$
(1)
Then, knowing $`P_2`$, $`A_j^2`$ and $`\phi _j^2`$ from the raw light curve we have removed the variability with a period $`P_2`$. For such prewhitened light curve we have computed the AoV spectrum which is displayed in Fig. 4b. The only frequency left is $`f_1=11.807`$ cycles/day corresponding to the period $`P_1=0.0847`$ days.

On the other hand having $`P_1`$, $`A_j^1`$ and $`\phi _j^1`$ one can make prewhitening of the light curve with period $`P_1`$. The result of such an operation is shown in Fig. 4c. In this case the only frequency that left is $`f_2=12.302`$ cycles/day which corresponds to the period $`P_1=0.0813`$ days.

Finally we have prewhitened our raw light curve of V1974 Cygni removing the variability with both periods $`P_1`$ and $`P_2`$. The AoV power spectrum of such prewhitened light curve is shown in Fig. 4d. There are no clear peaks in this graph indicating that $`f_1`$ and $`f_2`$ are sole frequencies present in the light curve of Nova Cygni 1992.

Knowing the values of both frequencies detected in the light curve of the star we were able to determine the amplitudes and shapes of light variations. The upper panel of Fig. 5 shows the raw light curve of V1974 Cygni phased with more distinct period $`P_1=0.08470`$ days. Due to the presence of the weaker signal characterized by period $`P_2=0.0813`$ days phased light curve shows significant scatter. Thus, based on equation (1), we removed the variations with $`P_2`$ and such a prewhitened light curve again phased with period $`P_1`$. The result is shown in the middle panel of Fig. 5. The points in this graph were additionally averaged in 0.025 phase bins. The solid line corresponds to the part of the equation (1) containing terms with $`P_1`$. The full amplitude of the modulation is $`0.090\pm 0.013`$ mag.

The lower panel of Fig. 5 shows a phased light curve with removed variations with period $`P_1`$ and its three harmonics. Again the points in this graph are the average values of the 0.025 phase bins. The full amplitude of this modulation is $`0.057\pm 0.015`$ mag.

## 5 O – C Diagrams

### 5.1 The 0.0813-day Period

The $`OC`$ analysis of times of extrema in a variable star is a very useful tool for tracing the changes of the period of the variable. Such an analysis, performed by DeYoung and Schmidt (1994), Semeniuk et al. (1995), Olech et al. (1996) and Skillman et al. (1997), showed that the 0.0813-day period of Nova Cygni 1992 is stable and can be interpreted as an orbital period of the binary system. Additionally Olech et al. (1996) pointed out that in 1994 this period underwent a slight change from 0.0812623 to 0.0812580 days, i.e. shortened by about 0.4 seconds. A possible mechanism of such a shortening of the orbital period may be the outflow of mass from the system through an outer Lagrangian point.

Now we can derive more definite conclusions using all maxima from the years 1993-1996 presented in papers by Semeniuk et al. (1995), Olech et al. (1996), Skillman et al. (1997) and from the years 1996-1997 presented in this paper.

Table 2 contains 69 times of maxima of the 0.0813-day modulation presented in the above-mentioned papers and in this work. The maximum from the night of October 4/5 1996 was determined as the average value of three peaks detected in our data. Because the amplitude of the modulation (see Fig. 2) was very high, we assumed that during this night the maxima of both periodicities were occurring at the same time, thus strengthening the modulations. Light variations detected in the light curve in July 1997 had a more complicated pattern, thus we decided to use it for only one maximum determination based on the $`P_2`$ term of the fit given by equation (1).

The 0.0813-day period behavior in the 1993-1997 time interval is illustrated in Fig. 6. The $`OC`$ residuals displayed in the figure and shown in Table 2 were calculated with the linear ephemeris:

$$\mathrm{HJD}_{\mathrm{max}}=2449535.4395+0.08125873E$$
(2)
$`\pm 0.0010\pm 0.00000023`$

Fig. 6 shows clearly that since the change in 1994 the orbital period of the binary system is stable and in the years 1994-1997 is equal to $`P_{orb}=0.08125873(23)`$ days.

### 5.2 The 0.085-day Period

Contrary to the 0.0813-day period, the longer period showed a more complicated behavior. Not present in the 1993 data, the 0.085-day period was the dominant peak in the light curve from the years 1994-1997 (Semeniuk et al. 1994, 1995, Olech et al. 1996, Skillman et al. 1997, Retter et al. 1997).

It was first demonstrated by Semeniuk et al. (1995) that the 0.085-day period is changing. Data collected in 1994 clearly showed a decrease of the period with a rate $`\dot{P}=7.8\times 10^7`$. In the beginning of 1995 the period started to increase (Olech et al. 1996) and was still increasing with $`\dot{P}=5\times 10^7`$ in 1996 (Skillman et al. 1997).

We decided to use all times of extrema of V1974 Cygni existing in the literature and determined in this work to trace the behavior of the 0.085-day period in detail. We used maxima given by Semeniuk et al (1994, 1995), Olech et. al (1996) and Retter et al. (1997). Skillman et al. (1997) tabulated only times of minima of 0.085-day periodicity but based on his Fig. 6 we transformed them into maxima adding the 0.45 phase shift to the minima.

Table 3 contains 139 times of maxima of the 0.085-day modulation presented in above-mentioned papers and determined in this work.

The 0.085-day period behavior in the 1994-1996 time interval is presented in Fig. 7. The $`OC`$ residuals displayed in the figure and shown in Table 3 were calculated with the linear ephemeris:

$$\mathrm{HJD}_{\mathrm{max}}=2449535.4320+0.085E$$
(3)
The fitted curve shown as a solid line corresponds to the ephemeris:

$$\mathrm{HJD}_{\mathrm{max}}=2449535.4406+0.08500527E1.435\times 10^8E^2+1.752\times 10^{12}E^3$$
(4)
Lack of the data between 1995 and 1996 seasons may suggest ambiguity in cycle count. In our opinion it is not in this case because equation (4) produces the value of the period equal to 122.66 min which is almost exactly the same as was reported by Skillman et al. (1997) from 1996 data i.e. $`P=122.67\pm 0.02`$ min. Cycle numerations in 1996 differing by plus or minus one from the scheme presented in Table 3 give the periods 122.63 and 122.74 min, respectively i.e. marginally consistent with data of Skillman et al. (1997).

Due to the lack of data in the almost 300-day interval between October 1996 and July 1997 a construction of the $`OC`$ diagram for the 1996-1997 data was impossible.

To trace the 0.085-day period behavior in 1997 we prewhitened the raw light curve of V1974 Cygni using second term of equation (1). From the resulting light curve we derived 8 times of minima of 0.085-day periodicity. They are shown in Table 4 together with $`OC`$ residuals computed from the following ephemeris:

$$\mathrm{HJD}_{\mathrm{max}}=2450634.3639+0.084632E$$
(5)
$`\pm 0.0068\pm 0.000080`$

Our 1997 observing run, which covers over 100 cycles, is too short to derive any conclusions about the 0.085-day period changes. Its mean value obtained from the $`OC`$ residuals is $`0.084632(80)`$ day i.e. $`121.87\pm 0.12`$ min which is consistent with $`122.0\pm 0.2`$ min value obtained from analysis of AoV power spectra.

As we mentioned earlier the 0.0813-day period was stable in the years 1994-1997 and its value was $`117.0126\pm 0.0003`$ min. The 0.085-day periodicity showed a more complicated behavior. It was decreasing form the middle of 1994 until the beginning of 1995. Semeniuk et al. (1995) and Skillman et al. (1997) reported a decrease as fast as $`\dot{P}=7.8\times 10^7`$ but our $`OC`$ diagram composed from all data available in the literature indicates a much slower decrease with a rate of $`\dot{P}=0.3\times 10^7`$.

From the beginning of 1995 the 0.085-day period started to increase and it was increasing until October 1996. Skillman et al. (1997) reported an increase with a rate $`\dot{P}=5\times 10^7`$ in 1996 but our $`OC`$ diagram shows that the rate was not higher than $`\dot{P}=1\times 10^7`$.

The most rapid change of the period was observed between October 1996 and July 1997. In this interval the period changed from 122.67 min (Skillman et. al. 1997) to 121.9 min. The $`\dot{P}`$ during this period of time could be as high as $`1.6\times 10^6`$ and certainly not slower than $`7\times 10^7`$.

## 6 Interpretation of the signals

Fig. 8 summarizes the results presented in this paper. It shows the evolution of both 0.085 and 0.081-day periodicities, their changes and evolution of the period of the disc precession assuming that V1974 Cygni is SU UMa type variable (see discussion below).

There is consensus in the literature concerning the interpretation of the 0.0813-day signal. Its discoverers, DeYoung and Schmidt (1994), suggested that it can be the orbital period of the binary system containing a white dwarf and a small main sequence star. The following papers by Semeniuk at al. (1995), Leibowitz et al. (1995), Olech et al. (1996), Skillman et al. (1997) and Retter et al. (1997) confirmed this conclusion.

Our investigation also confirmed the previous results. As was shown by Semeniuk et al. (1995) the value of the orbital period in 1993 was 0.0812623 days. In the first part of 1994 the period decreased by 0.4 seconds as was interpreted as a possible mass loss from the binary system. We have shown that in the years 1994-1997 the 0.085-day period was stable and its exact value was equal to $`0.08125873(23)`$ days ($`117.0126\pm 0.0003`$ min) supporting the hypothesis that it is indeed the orbital period of the binary system.

On the other hand there are two possible interpretations of the longer 0.085-day period. Semeniuk et al. (1994) who discovered it, suggested that it can be the rotation period of the white dwarf. According to that paper the evolution of the V1974 Cygni resembled the behavior of intermediate polars and particularly the famous Nova V1500 Cygni. Prior to the nova explosion, V1974 Cygni was a magnetized, synchronously rotating AM Her type binary. During outburst the spin/orbit synchronism was broken and the 0.085-day period, decreasing in 1994, may be the spin period of the white dwarf now evolving into synchronism with the orbital cycle. According to Semeniuk et al. (1995) the predicted time scale for equalization of both periods was $`13\pm 1`$ yr.

Observations from 1995 analyzed by Olech et al. (1996) showed, however, that the 0.085-day period stopped to decrease and since the middle of July 1995 seems to be increasing. Olech et al. (1996) pointed out that changes of the time derivative of spin periods are known in intermediate polars. FO Aqr is a good example of such objects (Osborne and Mukai 1989, Kruszewski and Semeniuk 1993). An additional argument for the spin nature of the 0.085-day period given by Olech et al. (1996) was a characteristic shape of the 1995 light curve of V1974 Cygni. A distinct feature of these light curves was a double structure of the maxima which is often observed in polars, possibly due to the existence of two magnetic poles (e.g. AM Her itself, Olson 1997)

Another explanation of the 0.085-day periodicity assumes that it might be the superhump period of a SU UMa type star (Leibowitz et al. 1995, Skillman et al. 1997, Retter et al. 1997). SU UMa stars are non-magnetized dwarf novae showing outbursts on a time scale of days and weeks and superoutbursts on a time scale of months. During the superoutburst a tidal instability causes the precession of the accretion disc. The presence of both precession and orbital periods causes the appearance of the third, so-called superhump period which is a few percent longer than the orbital period. In this interpretation V1974 Cygni is not a typical SU UMa star but rather a member of a newly discovered subgroup of these variables called permanent superhumpers (Skillman and Patterson 1993). These systems are characterized by mass transfer rates as high as $`\dot{M}=6\times 10^{16}`$ g s-1 (Osaki 1996) and thus are understood as nova-like stars below the period gap that are thermally stable but tidally unstable, i.e. being in permanent superoutburst.

Assuming that V1974 Cyg is a permanent superhumper we can directly compute its period of the disc precession which is a beat period between orbital and superhump periods. The behavior of the precession period in the years 1994-1997 is shown in lower panel of Fig. 8. Such rapid changes of the accretion disc are consistent with recent models of SU UMa systems (see for example Kornet and Różyczka 2000 or Buat-Menart and Hameury 2002)

About ten years ago the diversity of behavior observed in cataclysmic variables seemed to be quite small. These variables were classified into few distinct groups which were clearly visible in the $`P_{orb}\dot{M}`$ diagram (see for example Fig. 3 of Osaki 1996). There was a group of SU UMa stars below the period gap and at longer orbital periods one can find three other groups of U Gem, Z Cam and nova-like stars, listed according to the increasing $`\dot{M}`$.

In the mid-1990s the situation started to become more complicated. The SU UMa group itself was divided into several subgroups as ER UMa stars (Kato and Kunjaya 1995) and WZ Sge stars (Pych and Olech 1995, Kato et al. 1996). Additionally Nogami et al. (1998a, 1998b) and Uemura et al. (2001) discovered the first dwarf novae in the middle of the period gap. The stars with short orbital periods and high accretion rates were found to be permanent superhumper systems (Skillman and Patterson 1993). Superhump-like variations were also detected in ultra-short orbital period systems containing helium secondary (e.g. AM CVn, V485 Cen, HP Lib or 1RXS J232953.9+062814 - Skillman et al. 1999, Olech 1997, Patterson et al. 2002, Skillman et al. 2002) and possibly a black hole primary (e.g. KV UMa, Uemura et al. 2002)

Permanent superhumps were also thought to be the reason for periodic light variations in other high accretion binaries on both sides of the period gap. At longer periods positive and negative superhumps were detected in nova-like variables (for example V592 Cas and V751 Cyg - Taylor et al. 1998, Patterson et al. 2001).

Finally, permanent superhumps were assumed as the main cause of short period variations of the light curves of classical nova systems lying on both sides of the periodic gap (e.g. V603 Aql, CP Pup and V1974 Cyg itself - Thomas 1993, Patterson and Warner 1998, Retter and Leibowitz 1998).

Retter and Leibowitz (1998) found that these three classical novae are characterized by high accretion rates ($`\dot{M}`$ from $`8\times 10^{16}`$ to $`8\times 10^{17}`$ g s-1) and thus their periodic variations could be interpreted as permanent superhumps. The most interesting case is V1974 Cygni itself. According to Retter and Leibowitz (1998), it should undergo transition from a permanent superhumper to the ordinary SU UMa star around the year 2000 (with an uncertainty of about 2 years).

To summarize, the latest results obtained for several cataclysmic variables strongly supports the hypothesis that 122-min periodicity observed in the light curve of V1974 Cygni is a superhump period. The period changes and the changes of the period derivative are known in both SU UMa stars and polars. But the temporal evolution of the 122-min period especially in 1996-1997 was quite rapid. Such rapid changes of the period (sometimes as high as $`\dot{P}`$ around $`10^5`$) arising from changes in the behavior of the accretion disc are observed in SU UMa systems. On the other hand the $`0.60.8`$ min change of the 122-min periodicity which occurred between 1996 and 1997 would be problematic from the point of view of the hypothesis assuming that we observe the rotation period of a magnetized white dwarf. Such a rapid change would require a change in the rotational kinetic energy of the white dwarf as large as $`10^{42}`$ ergs. This argument was also given by Skillman et al. (1997) based on the period changes of V1974 Cygni in 1995-1996. In our work, however, we showed that these changes were not as abrupt as was previously thought, thus the change of the rotational energy was not so large. In the end of 1996 the 122-min period started to decrease more rapidly and thus we again should expect a change of the rotational energy at the level of $`10^{42}`$ ergs which is higher than the entire energy radiated by the star in the October 1996 - July 1997 interval.

Finally, we conclude that the explanation of the 122-min periodicity involving the permanent superhumps phenomenon is more consistent with observations than assumption that we are observing the rotation of a white dwarf. However, the beginning of the new millennium is the ideal time for checking both hypotheses. The lack of the 122-min periodicity in the present light curve of Nova Cygni 1992 would support the superhump explanation indicating that the object became an ordinary SU UMa star and the presence of this period in circular polarization would validate the connection of the period with white dwarf rotation.

Acknowledgments. We would like to thank to Prof. Józef Smak and Chris O’Connor for reading and commenting on the manuscript. Author was supported by the Polish State Committe for Scientific Research (KBN), grant 2P03D02422, and a Foundation for Polish Science stipend for young scientists.


text:
# Schrödinger equation for current carrying states

## I Introduction

The remarkable miniaturization of conventional microelectronic devices has been witnessed over the last decades. If the miniaturization trend is to continue, elements of microelectronic circuits, e.g. transistors and contacts, will soon shrink to single molecules or even atoms. In the past few years, there have been considerable advances towards fabrications and experimental studies of mono-molecular electronic devices joachim2000 . This activity has been largely spurred on by development of experimental techniques to form atomic scale electrical contacts such as scanning tunneling microscopy mirkin92  and mechanically controllable break junctions ruiten96 . Measurements of current-voltage (I-V) characteristics have been recently performed on a single benzene-1,4-dithiolate molecule reed97 , buckminsterfullerene joachim97 , individual atoms scheer98 , nanoclusters wildoer98  and carbon nanotubes martel98 . These initial observations of conduction in molecules were followed quickly by experimental demonstration of the simplest devices: molecular diodes metzger2000 , switches chen99  and memory elements reed2001  although molecular equivalent of a transistor is yet to be discovered.

Several theoretical approaches have been recently developed to describe electronic conductivity of molecules nitzan2001 . One of the widely used methods is based on the combination of extended Hückel model or density functional theory to treat the molecular electronic structure, an electronic reservoir model for the contacts and linear response approximation to compute the electronic current as a function of the applied voltage mujica94 ; tian98 ; hall2000 ; gutierrez2001 ; derosa2001 . The most recent applications used Lippmann-Schwinger or non-equilibrium Green’s functions formalisms plus “jellium” model for the contacts and for the molecular-contact interactions diventra2000 ; taylor2001 . In all these approaches, the system is partitioned into three parts: molecular wire and two electron reservoirs to represent metallic contacts. The use of the model dependent reservoir-molecule interaction as one of the most pivotal ingredient diminishes the ability of these methods and prevents from a development of predictive and quantitative techniques for the electronic transport problem. Recent density functional theory based calculation demonstrates that a molecular wire is bonded to a gold contact by the thiol bridge which is known to be the strong covalent interaction johansson2000 . The empirically treated coupling between metal contacts and the molecular wire, which serves as the interaction for the linear response calculations and represents the kernel of the corresponding Lippmann-Schwinger equation, is at least strong as intermolecular interaction. Recent papers also indicate the exact nature of the molecule-contact bonding is critical in predicting the correct order of magnitude for the current with applied bias diventra2000 .

The miniaturization shrinks the size of the system to a level when the problem can be addressed by ab initio electronic structure methods. As an example consider a standard prototype molecular device system with a single benzene-1,4-dithiolate molecule covalently bonded to two gold electrodes reed97 . The problem can be reduced to the system where the molecule is bonded to two small Au clusters derosa2001 . A gold cluster with reasonably small edge influence on the molecular-surface bonding contains 10 atoms johansson2000 . With a single benzene-1,4-dithiolate molecule as the molecular wire the contact-wire-contact system consists of just 32 atoms. Systems of such size can be routinely simulated by standard quantum chemistry codes. The main obstacle for the application of the first principle electronic structure methods is not the size of the system but modeling of the electronic current and interpretation of the applied voltage. Suppose that we have solved Hartree-Fock or Kohn-Sham equations for a molecular device. The ground state wave function is real in absence of a magnetic field and therefore the current density computed from this wave function is zero everywhere in the sample. Another important issue is an interpretation of the applied voltage. Virtually all present theoretical approaches are pivotally relied on simplified, noninteracting electrons type models for the metal contacts where the external voltage bias can be straightforwardly associated with the Fermi energies shift between left and right electronic baths. One of the reasons for use electronic reservoir models to represent the metallic contacts is principle theoretical difficulties related to a rigorous quantum mechanical definition of the applied voltage.

In this paper, we advocate an alternative approach where the use of the electronic reservoir to model contacts is completely avoided from the very beginning and a quantum system with current is treated effectively as a closed  system. We begin by describing the variational principle for a quantum system with current. We then derive the Schrödinger equation with probability flux which provides exact wavefunctions on the manifold of a given current density. We next specify two-terminal molecular device steady current subsidiary condition and develop the Schrödinger equation for a two-terminal molecular device. We discuss the use of the energy ballance condition to compute current-voltage characteristics of quantum systems. We then derive the Landauer formula for non-interacting charge carriers within the formalism. Finally, we demonstrate the salient features of the approach through the numerical solution of the self-consistent Schrodinger equation for one-dimensional harmonic oscillator with constrained current.

## II Schrödinger equation with probability current

### II.1 Variational principle

We start with the nonrelativistic Hamiltonian for a particle in the external scalar potential $`v(𝐫)`$ (atomic units are used throughout the paper unless otherwise mentioned):

$$H=\frac{1}{2}^2+v(𝐫).$$
(1)
We define density and paramagnetic current density as:

$`\rho (𝐫)`$ $`=`$ $`\psi (𝐫)\psi ^{}(𝐫)`$ (2)
$`𝐣(𝐫)`$ $`=`$ $`\text{Im}\left\{\psi ^{}(𝐫)\psi (𝐫)\right\}`$ (3)
where $`\psi (𝐫)`$ is not generally an eigenvector of the Hamiltonian $`H`$. We will specify an eigenproblem for the $`\psi (𝐫)`$ in the next section. The physical current is equivalent to the paramagnetic current densities in absence of an external magnetic field.

An eigenvector for a bound state of a real Hamiltonian can be always made real by a change of the phase of the wave function. The Hamiltonian $`H`$ is real therefore if we just solve standard Schrödinger equation with the Hamiltonian $`H`$, we will obtain real wave function $`\psi (𝐫)`$ and the current density $`𝐣(𝐫)`$ will be automatically zero. The standard way to get a complex wavefunction out of a real Hamiltonian is to assume certain boundary conditions. These boundary conditions could be an assumption of an incoming plane wave in scattering theory or periodic boundary conditions for a translationally invariant system. To solve a scattering problem for a molecular device is possible only at the cost of major unphysical approximations, e.g. via considering the molecular device as a singe scattering center. Another approach is to assume periodic boundary conditions. The periodic boundary conditions prevent one from studying any situations where the change of chemical potential across the system is of finite magnitude because the chemical potential must be also periodic in space. The local chemical potential is the decreasing function in the direction of the current and therefore the system can not be treated as translationally invariant at least along the current axis.

We have developed an alternative approach where a current carrying quantum system is treated effectively as a closed  system with a subsidiary condition for the generally discontinuous current density distribution. We start our derivation with the assumption that current carrying states of an open quantum system can be described in terms of wave functions of an effectively closed quantum system if the current is maintained via explicit constraint. This assumption is a simplified, particular case of the standard assumption in statistical mechanics: a system coupled to a thermal bath can be described with statistical ensembles of a closed system in which the system-bath interaction is not explicitely considered, but enters only via Lagrange multipliers and subsidiary constraints. It has been recently demonstrated that a steady current open mesoscopic system can be described in terms of an ensemble of carrier states in the closed mesoscopic device itself if the steady current is additionally constrained heinonen93 ; johnson95 .

We require that the current density is constrained to be specified function $`𝐈(𝐫)`$. It results in the subsidiary constraint equation:

$$𝐣(𝐫)𝐈(𝐫)=0.$$
(4)
The expression (4) generally can not be obtained by differentiating of a functional of $`\psi `$ and $`\psi ^{}`$ with respect to $`𝐫`$. It means that the constraint eq.(4) is nonholonomic courand . The variational problem is to minimize the total energy with respect to $`\psi `$ subject to the nonholonomic constraint for the current density (4). Two types of treatment are possible: we can either choose $`\psi (𝐫)`$ in which the constraint is implicit or we can impose the constraint explicitely by using the Lagrange multiplier. We will follow the second approach which was recently suggested by Kosov and Greer kosov2001 , elimination the part of the wavefunction by imposing the explicit constraint. Although the constraint formulated as eq.(4) is vector and nonholonomic it can be included into a variational functional via pointwise, vector Lagrange multipliers $`𝐚(𝐫)`$ kosov2001 .

Within the Lagrange multiplier approach the constraint is explicitly achieved via introduction of the auxiliary functional:

$$\mathrm{\Omega }[\psi ]=\psi |H|\psi E\left(\psi |\psi 1\right)𝐚(𝐫)\left(𝐣(𝐫)𝐈(𝐫)\right)𝑑𝐫.$$
(5)
The first two terms are standard in the variational derivation of the Schrödinger equation with the first giving the total energy and the second is introduced to maintain the normalization of the wave function. The third term with the vector Lagrange multipliers $`𝐚(𝐫)`$ has been introduced to impose the subsidiary constraint for the current density.

A use of the variational principle can be justified only for closed quantum systems. The choice of $`𝐈(𝐫)`$ can not be obviously completely arbitrary to maintain the closed system boundary condition. The localized in space quantum system will be considered closed if, first, there are no in- and out- probability flows in the system; and, second, if the wave function $`\psi `$ is an eigenvector of a Hermitian operator (it prevents the complex eigenenergies and therefore effectively closes the system). We will demonstrate in the next section that both these conditions are satisfied within our approach.

### II.2 Self-consistent Schrödinger equation

The minimization of the total energy subject to the subsidiary constraint for the current density (4), i.e the variation of the auxiliary functional $`\mathrm{\Omega }[\psi ]`$, results into the following non-linear Schrödinger equation:

$$\left(\frac{1}{2}\mathrm{\Delta }+v(𝐫)\frac{1}{2i}[𝐚(𝐫),]_+\right)\psi (𝐫)=E\psi (𝐫).$$
(6)
The anti-commutator term,

$$\frac{1}{2i}[,𝐚(𝐫)]_+=\frac{1}{2i}\left(𝐚(𝐫)+2𝐚(𝐫)\right),$$
(7)
is the additional completely imaginary potential arising directly from the constraint on the current density (consult with appendix A for the corresponding functional derivatives ). This additional constraint potential forces the wave function $`\psi (𝐫)`$ to be irremovably complex and insures that $`\psi (𝐫)`$ produces the required current density $`𝐈(𝐫)`$.

A physical interpretation of the Lagrange multiplier $`𝐚(𝐫)`$ can be made more explicit if we re-write eq.(6) in the following form:

$$\left[\frac{1}{2}\left(i𝐚(𝐫)\right)^2+\frac{1}{2}𝐚(𝐫)𝐚(𝐫)+v(𝐫)\right]\psi (𝐫)=E\psi (𝐫).$$
(8)
It can be now established that the current constraint has introduced the terms which are mathematically equivalent to an external vector potential. But unlike an external vector potential, the Lagrange multiplier $`𝐚(𝐫)`$ does depend upon the wave function $`\psi (𝐫)`$ and eq.(8) becomes nonlinear and must be solved self-consistently. The eigenvalue problem eq.(8) is generated by Hermitian operator since the Lagrange multiplier $`𝐚(𝐫)`$ is real. Being a solution of the Hermitian eigenproblem the eigenenergy $`E`$ is real and the wavefunctions $`\psi `$ form a complete set in Hilbert space. Within our approach we deal with the localized in space system with the wave functions generated by the hermitian eigenproblem (8). Such system is by definition closed and therefore the use of the variational principle is formally justified.

The Schrödinger equation with current (8) is not yet in a form allowing for a solution to be found as the Lagrange multiplier $`𝐚(𝐫)`$ is not known yet. We next seek a solution for the Lagrangian multiplier $`𝐚(𝐫)`$. In the developing of this method we have found out that the most significant source of numerical instabilities in the solution of the nonlinear Schrödinger equation (8) arises from the determination of the Lagrange multiplier $`𝐚(𝐫)`$. The described below procedure obtaining of the Lagrange multiplier has been found to be the most robust. Multiplying eq. (6) from the left on the $`\psi (𝐫)^{}`$ and subtracting from the obtained equation its own Hermitian conjugate, yields

$$𝐣(𝐫)=(𝐚(𝐫)\rho (𝐫)).$$
(9)
Equation (9) can be resolved for $`𝐚(𝐫)`$ with the assumption that both the density $`\rho (𝐫)`$ and the current density $`𝐣(𝐫)`$ decay to zero at infinity. In this case the direct integration of eq.(9) results into the following expression for the Lagrange multiplier $`𝐚(𝐫)`$:

$$𝐚(𝐫)=\frac{𝐣(𝐫)}{\rho (𝐫)}$$
(10)
Substitution of the constraint for the current density eq.(4) into eq.(10) yields the simple expression in which the Lagrange multiplier depends upon the wave function only via the density in the denominator:

$$𝐚(𝐫)=\frac{𝐈(𝐫)}{\rho (𝐫)}$$
(11)
There is a direct similarity between equation for the Lagrange multiplier (10) and the definition of the velocity of a quantum trajectory within de Broglie-Bohm causal description of quantum mechanics holland . Within the de Broglie-Bohm dynamics the velocity $`𝐮(𝐫,t)`$ of a trajectory at a given point is computed as bittner2000 :

$$𝐮(𝐫,t)=\frac{𝐣(𝐫,t)}{\rho (𝐫,t)},$$
(12)
therefore, at least for stationary systems the vector Lagrange multiplier $`𝐚(𝐫)`$ can be directly related to the velocity for a corresponding de Broglie-Bohm trajectory. By this analogy we can make the physical interpretation of the Lagrange multiplier more explicit. In the de Broglie-Bohm quantum dynamics to maintain the constant flux at the given density profile one needs to adjust the local trajectory velocities. The vector Lagrange multiplier appears to play the simular role in the Schrodinger equation with constrained current. Likewise for a point with the non-zero trajectory velocity, the current density has finite value at a point with the non-zero vector Lagrange multiplier. And likewise the de Broglie-Bohm velocity, a product of the vector Lagrange multiplier on the probability density yields the current density.

With this choice of the Lagrange multiplier (11) we obtain the Schrödinger equation with the current $`𝐈(𝐫)`$ in the self-consistent and closed form:

$$\left[\frac{1}{2}\left(i\frac{𝐈(𝐫)}{\rho (𝐫)}\right)^2+\frac{1}{2}\frac{I^2(𝐫)}{\rho ^2(𝐫)}+v(𝐫)\right]\psi (𝐫)=E\psi (𝐫).$$
(13)
The derivation of the self-consistent Schrödinger equation eq.(13) composes one of the central results of the paper.

## III Schrödinger equation for a two terminal molecular device

The desire to develop the basic Schrödinger equation for calculations of $`IV`$ characteristics of two terminal molecular devices is the impetus for the present study. To achieve this aim, we specify the molecular wire constraint on the current density:

$$𝑑y𝑑zj_x(𝐫)=I_x\theta (xL,Rx),$$
(14)
where $`I_x`$ is the steady current through the molecular device. The step function,

$`\theta (xL,Rx)=\theta (xL)\theta (Rx),`$ (15)
equals to 1 if $`LxR`$ and zero otherwise. $`L`$ and $`R`$ are the left and right boundaries for the molecular device. Within this description, net current flow is aligned along the $`x`$axis and there are no in- and out flow: the current starts at the left boundary L and completely absorbed at the right R. With this constraint it is not required to fix the full vector $`𝐣(𝐫)`$ everywhere in the sample. Only the net current flow across a cross section $`𝑑y𝑑zj_x(𝐫)`$ is constrained, and this quantity is readily available experimentally. This is a simple geometric arrangement to specify current flow for the two terminal molecular device and can be extended to more complex device topologies, e.g. for the three terminal devices where three tips are connected to a single molecule diventra2000\_b .

Likewise, as for the derivation of the exact Schrödinger equation in the preceding section, with this choice of the constraint (14) the minimization of the auxiliary functional (5) now yields the Schrödinger equation for the steady current two terminal molecular device:

$$\left(\frac{1}{2}\mathrm{\Delta }+v(𝐫)\frac{1}{2i}[a_x(x),\frac{}{x}]_+\right)\psi (𝐫)=E\psi (𝐫).$$
(16)
The expression for the Lagrange multiplier given by eq.(11) can be straightforwardly extended to the steady current molecular device, with the result is

$$a_x(x)=\frac{I_x}{\rho _{yz}(x)}\theta (xL,Rx),$$
(17)
where we have introduced the quantity

$$\rho _{yz}(x)=𝑑y𝑑z\rho (𝐫).$$
(18)
With the expression for the Lagrange multiplier $`a_x(x)`$ in hands, we can complete the Schrödinger equation for the two-terminal molecular device with steady current $`I_x(x)`$

$$\left(\frac{1}{2}\mathrm{\Delta }+v(𝐫)\frac{I_x}{2i}[\frac{\theta (xL,Rx)}{\rho _{yz}(x)},\frac{}{x}]_+\right)\psi (𝐫)=E\psi (𝐫).$$
(19)
Some aspects of eq.(19) deserve special discussion. The imaginary potential which enforces the wave function to produce steady current $`I_x`$ depends upon the $`x`$ coordinate only. The range of this potential is restricted by the molecular device boundaries. One of the terms in this potential is proportional to the derivative of the step function,

$$\frac{}{x}\theta (xL,Rx)=\delta (xL)\delta (xR),$$
(20)
and therefore yields the singular imaginary $`\delta `$-function potential on the device boundaries. It is known from the standard quantum mechanics textbooks fluege71  that a $`\delta `$-function potential results in the discontinuity of the first derivative of the corresponding solutions of the Schrödinger equation. Likewise, the imaginary component of the eigenfunction $`\psi (𝐫)`$ of the Schrödinger equation for a molecular device (19) has the discontinuous first derivative at the points $`x=L`$ and $`x=R`$.

The final comment regarding the computation of the I-V characteristic of a two-terminal molecular device is in due order. The current carrying states obtained as the solutions of the molecular device Schrödinger equation (19) must now be associated with current-voltage characteristics available experimentally. The rigorous definition of the applied voltage is a controversial issue if one does not invoke to noninteracting electron reservoirs to represent contacts. We propose the following macroscopic definition of the applied voltage which has the direct relation with the quantities measured experimentally. The energy increase from the applied voltage $`U`$ is computed by the numerical integration with an assumption of linear voltage drop between the left and right boundaries:

$$\mathrm{\Delta }E_d(U)=\frac{U}{RL}\underset{L}{\overset{R}{}}𝑑x(Rx)\rho _{yz}(x).$$
(21)
We compute the energy of the non-active, i.e. with zero current, molecular device and compute the total energy of the same molecular device at an experimentally given current $`I_x`$. Then the energy increase for the establishment of the current carrying state, i.e. the energy difference

$$\mathrm{\Delta }E_d(I_x)=E_d(I_x)E_d(I_x=0)$$
(22)
is calculated. Because of the energy ballance condition the $`\mathrm{\Delta }E_d(U)`$ must be equal to the $`\mathrm{\Delta }E_d(I_x)`$ and it results into the following equation for the applied voltage bias:

$$U=\frac{(RL)}{\underset{L}{\overset{R}{}}𝑑x(Rx)\rho _{yz}(x)}\mathrm{\Delta }E_d(I_x).$$
(23)
The eigenenergy $`E`$ of the Schrödinger equation (19) should not be confused with the energy of molecular device. The eigenenergy $`E`$ serves only as the Lagrange multiplier to enforce orthogonality of the wave function. The energy of the molecular device can be obtained as an expectation value of the real physical Hamiltonian $`H`$ with the integration restricted by the molecular device region:

$$E_d(I_x)=𝑑y𝑑z\underset{L}{\overset{R}{}}𝑑x\psi ^{}(𝐫)\left(\frac{1}{2}\mathrm{\Delta }+v(𝐫)\right)\psi (𝐫),$$
(24)
where $`\psi (𝐫)`$ is the solution of the Schrödinger equation (19) with the given $`I_x`$. Since the wavefunction $`\psi (𝐫)`$ is complex and the integration (24) is not performed over whole space, it is not obvious from the expression (24) that the device energy $`E_d(I_x)`$ must be real. We will show now that the device energy $`E_d(I_x)`$ is real if the wavefunction $`\psi (𝐫)`$ produces the constrained current density (14). To prove this we subtract from eq.(24) its own complex conjugate and after a bit of algebra arrive to the following expression:

$$E_d(I_x)E_d(I_x)^{}=i𝑑y𝑑z\underset{L}{\overset{R}{}}𝑑x𝐣(𝐫)=i_S𝑑S𝐧𝐣(𝐫).$$
(25)
The length of the quantum wire is aligned along the $`x`$axis and the surface integral does the integration over the surface of the rectangular box with one $`yz`$ face at $`x=L`$, another $`yz`$ face at $`x=R`$ and any of the $`xy`$ or $`xz`$ faces are removed at the infinity. Within our constraint, there is a positive current flowing into the $`yz`$ face at $`x=L`$ and exactly the same positive current flowing out of the $`yz`$ face at $`x=R`$. There is no current flow out of the region in any of the $`xy`$ or $`xz`$ faces. Therefore the surface integral vanishes

$$_S𝑑S𝐧𝐣(𝐫)=0$$
(26)
and the device energy $`E_d(I_x)`$ defined via eq.(24) becomes real.

## IV Example calculations

### IV.1 Conductance of an ideal one-dimensional lead

In order to establish connection with standard approaches and test the method we derive the Landauer formula landauer70  using the Schrödinger equation for current carrying states eq.(13) as the starting point. Consider noninteracting electrons moving along a one-dimensional wire without any scatterer. The electronic density does not depend upon $`x`$ yielding the following Schrödinger equation with steady current $`I`$:

$$\left[\frac{1}{2}\frac{d^2}{dx^2}+i\frac{I}{\rho }\frac{d}{dx}\right]\psi _k(x)=E_k\psi _k(x).$$
(27)
Equation (27) has the plane-wave solution

$$\psi _k(x)\mathrm{exp}(ikx)$$
(28)
with eigenenergy

$$E_k=\frac{k^2}{2}\frac{Ik}{\rho }.$$
(29)
The similar dispersion relation as given by eq.(29) has been obtained within the same approach for uniform electron gas with applied direct current kosov2001 . For the small voltage and temperature transmissions of electrons take place only through the states in close vicinity of the Fermi level. There are two vectors k, one is positive and one is negative, which correspond to the Fermi energy:

$`k_+`$ $`=`$ $`{\displaystyle \frac{I}{\rho }}+\sqrt{2E_F+\left({\displaystyle \frac{I}{\rho }}\right)^2},`$ (30)
$`k_{}`$ $`=`$ $`{\displaystyle \frac{I}{\rho }}\sqrt{2E_F+\left({\displaystyle \frac{I}{\rho }}\right)^2}.`$ (31)
Following the Landauer approach we associate the voltage drop (at the zero temperature) with the gap between the singe-electron energies of the electrons moving in the direction of the current ($`k>0`$) and electrons moving in the opposite direction ($`k<0`$):

$$U=\frac{k_+^2}{2}\frac{k_{}^2}{2}.$$
(32)
Substituting the expressions for the $`k_+`$ and $`k_{}`$ eqs.(30,31) into the Landauer’s definition of the applied voltage (32), we readily find

$$I=GU,$$
(33)
with the conductance

$$G=\frac{\rho }{2\sqrt{2E_F+\left(\frac{I}{\rho }\right)^2}}.$$
(34)
Then, the use of the standard relation between the Fermi momentum $`k_F`$ and the electronic density $`\rho =k_F/\pi `$ along with the small current expansion of the square root eq.(34) results into the Landauer formula for a single transport channel of an ideal one-dimensional lead datta95 :

$$G=\frac{1}{2\pi }.$$
(35)
If we take into account the spin degeneracy of real electrons, the conductance becomes multiplied by 2. By accounting the degeneracy and converting the conductance eq.(35) to SI units the standard value $`G=e^2/(\pi \mathrm{})=(12.9K\mathrm{\Omega })^1`$ is obtained. Re-deriving the standard result we demonstrate that in the low current and zero temperature limits our method is equivalent to the traditional approaches which are based on the Landauer formula.

### IV.2 Numerical solution of the Schrödinger equation with current

We demonstrate applications of the method with numerical examples. The numerical solution of a self-consistent Schrödinger equation is often performed within a finite basis set expansion for the wave function with subsequent diagonalization of the Hamiltonian matrix until the self consistency is reached. We have found out that it is computational convenient to form the complete, real basis set from the solutions of the zero-current Schrödinger equation:

$$\left(\frac{1}{2}\mathrm{\Delta }+v(𝐫)\right)\varphi _\mu (𝐫)=\epsilon _\mu \varphi _\mu (𝐫).$$
(36)
Then we expand the wave function of a current currying state on this basis

$$\psi (𝐫)=\underset{\mu }{}C_\mu \varphi _\mu (𝐫),$$
(37)
the coefficients $`C_\mu `$ must be complex for the wave function which produces a non-zero current. With the expansion (37) we arrive to the following complex, Hermitian eigenproblem:

$$\underset{\mu }{}\left(\epsilon _\mu \delta _{\mu \nu }iw_{\nu \mu }(𝐈)\right)C_\mu =EC_\nu ,$$
(38)
where the current dependent matrix element $`w_{\nu \mu }(𝐈)`$ is computed in appendix B and is given by the following integral:

$$w_{\nu \mu }(𝐈)=\frac{1}{2}𝑑𝐫\frac{𝐈(𝐫)}{\rho (𝐫)}\left(\varphi _m(𝐫)\varphi _n^{}(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)\right).$$
(39)
The matrix elements $`w_{\nu \mu }(𝐈)`$ are real and depend on the expansion coefficients $`C_\mu `$ via the density $`\rho (𝐫)`$ thereby making the eigenvalue problem (38) nonlinear. The eigenproblem (38) is Hermitian because the matrix element $`w_{\nu \mu }(𝐈)`$ changes its sign when we swap $`\mu `$ and $`\nu `$ indeces:

$$w_{\nu \mu }(𝐈)=w_{\mu \nu }(𝐈).$$
(40)
The numerical calculations were carried out for the one-dimensional system. With an eye on molecular device simulations we solve the self-consistent Schrödinger equation for the step function constraint as a primary example:

$$j(x)=I\theta (xL,Rx)$$
(41)
The external potential was taken in the harmonic oscillator form and we put $`L=2.0`$ and $`R=2.0`$ as left and right boundaries for the current carrying part of the system in all our calculations. Our basis set was built up from the $`25`$ eigenvectors of the zero current Hamiltonian. The self-consistency was assumed to be achieved when the root mean square deviation of the wave function vectors of the two subsequent iterations was less then $`10^8`$.

First, we plot the Lagrange multiplier as a function of the coordinate on the figure 1. The Lagrange multiplier is discontinuous at the system boundaries. The Lagrange multiplier is vanishing out of the system region and it is the inverse function of the density $`\rho (x)`$ inside the system. As it can be expected from the analogy with the de Broglie-Bohm quantum hydrodynamics (12) (the lower the density, the larger the local velocity) the Lagrange multiplier reaches its maximum values in the regions of depleted density and has the minimal value at the point of maximal density.

In our next example we studied the deviation of the current carrying wave function $`\psi (𝐫)`$ ($`I=0.1`$) from the zero-current ground state wave function $`\varphi _0(𝐫)`$. The real component is depicted on the upper panel of the fig.2. The shape and magnitude of the real component do not undergo significant transformation with the establishment of the current carrying state. Although the admixture of the higher terms in the expansion (37) results in the small damping oscillations on the tails of the real component of the wave function. The imaginary part of the wave function ( plotted on the lower panel of the fig.2) is identically zero when the current is vanished. The imaginary component reaches its negative minimum at the left boundary $`L=2.0`$ of the current carrying part of the system and then it discontinuously changes the sign of the first derivative and undergoes almost linear increase until it reaches its maximum at the right boundary $`R=2.0`$. When the coordinate $`x`$ goes out of the non-zero current region the imaginary part of the wave function rapidly decays to zero value.

The establishment of the imaginary part of the wave function directly indicates the development of the current carrying state. We studied transformations of the imaginary component as a function of the increased current. The imaginary components of the wave functions are plotted on the fig. 3 for three representative steady currents: $`I=0.01`$, $`I=0.5`$ and $`I=0.1`$. At the very small current $`I=0.01`$ the imaginary wave function is almost zero everywhere except the only small negative minimum at the left boundary and with the positive maximum at the right. The magnitudes of the imaginary component at the extremums are nearly linearly proportional to the value of the constrained current. The slope of the imaginary component, i.e. the derivative, is always positive and almost constant between the left and right boundaries of the system. The imaginary part becomes zero when the real component reaches its maximum, i.e. at the point $`x=0.0`$. Although the total wave function is real at the point $`x=0.0`$ the current is not vanished because the derivative of the wave function remains complex.

Now we turn to the calculations of a “device” energy required to establish the current carrying state. Being the Lagrangian multiplier to maintain the wavefunction orthonormality, the eigenenergy $`E`$ of the eigenproblem (38) is not directly related to the energy of the system. The energy increase due to the establishment of the current carrying state is defined via the space restricted expectation value of the Hamiltonian (24):

$$E_d(I)=\underset{L}{\overset{R}{}}𝑑x\psi ^{}(x)\left(\frac{1}{2}\frac{d^2}{dx^2}+v(x)\right)\psi (x)=\underset{\mu \nu }{}C_\mu ^{}C_\nu \epsilon _\nu \underset{L}{\overset{R}{}}𝑑x\varphi _\mu (x)\varphi _\nu (x).$$
(42)
The table 1 contains the energy $`E_d(I)`$ and the energy difference between the nonactive, i.e. zero-current, and current $`I`$ carrying state as a function of the increased current. The “device” energy is gradually increased with increasing current. As it can be deduced from the table 1, we need to pump into the system $`0.14a.u.`$ of energy to establish the current $`I=0.1a.u.`$ carrying state. A simple fit of the data from the table 1 provides only the quadratic and cubic dependences of the “device” energy upon the steady current (for $`I0.1`$):

$$E_d(I)=E_d(I=0)+17.4I^236.5I^3,$$
(43)
thereby demonstrating that the term linear in the current has the minor importance for the energetics of the current carrying quantum system and the main energy contribution comes from the $`I^2`$ term. This energy dependence is generally agreed with the current-density functional theory of inhomogeneous interacting electron gas where the current dependent correction to the standard exchange-correlation energy functional is also mainly proportional to $`I^2`$ vignale87 .

## V Conclusions

In this paper, we have given a variational formulation of the Schrödinger equation with non-zero current. The Schrödinger equation with current is derived via a constrained minimization of the total energy with a subsidiary condition for the current density. The subsidiary condition for the current density is maintained during the course of variation by a vector, pointwise Lagrange multiplier. An explicit elimination of this Lagrange multiplier results into the closed, self-consistent form of Schrödinger equation with current. We showed that the current carrying states are eigenvectors of the complex, Hermitian operator. The formulation has been developed for general current density topologies and then specified for the case of a two-terminal molecular device. We showed how the energy ballance condition can be used for rigorous computation of I-V characteristics of molecular devices: the energy increase due to the establishment of a current carrying state is associated with the applied voltage bias. We demonstrate the salient features of the approach by re-deriving the Landauer formula for the conductance of an ideal one-dimensional lead and through a solution of the one-dimensional Schrödinger equation with fixed current density. The complex, Hermitian matrix was formed on a basis of the eigenvectors of the zero-current Hamiltonian. The numerical solution was achieved via the self-consistent solution of the complex, Hermitian eigenproblem. We found out that a basis set formed from the zero current Hamiltonian eigenvectors provides a robust algorithm for the self-consistent convergence. Our method has the great compatability with standard electronic structure methods which are also based on a variational principle e.g. Hartree-Fock or configuration interaction approximations.

###### Acknowledgements.

The author wishes to thank G. Stock for stimulating discussions, critical reading the manuscript and his continuous support of the project. The author also thanks E.R.Bittner, J.C.Greer and A.Nitzan for valuable discussions.

## Appendix A Functional derivatives of the constraint functional

In this Appendix we compute the functional derivative of the constraint functional. The constraint functional has the following form

$$\mathrm{\Lambda }[\psi ]=𝑑𝐫𝐚(𝐫)(𝐣(𝐫)𝐈(𝐫)).$$
(44)
The direct variation of the constraint functional with respect to the $`\psi ^{}(𝐫)`$ yields:

$`{\displaystyle \frac{\delta \mathrm{\Lambda }[\psi ]}{\delta \psi ^{}(𝐫)}}`$ $`=`$ $`{\displaystyle 𝑑𝐫^{}\frac{\delta \mathrm{\Lambda }}{\delta 𝐣(𝐫^{})}\frac{\delta 𝐣(𝐫^{})}{\delta \psi ^{}(𝐫)}}={\displaystyle 𝑑𝐫^{}𝐚(𝐫^{})\frac{\delta 𝐣(𝐫^{})}{\delta \psi ^{}(𝐫)}}`$
$`=`$ $`{\displaystyle \frac{1}{2i}}{\displaystyle 𝑑𝐫^{}𝐚(𝐫^{})\frac{\delta }{\delta \psi ^{}(𝐫)}\left\{\psi ^{}(𝐫^{})\psi (𝐫^{})\psi (𝐫^{})\psi ^{}(𝐫^{})\right\}}`$
$`=`$ $`{\displaystyle \frac{1}{2i}}𝐚(𝐫)\psi (𝐫){\displaystyle \frac{1}{2i}}{\displaystyle \frac{\delta }{\delta \psi ^{}(𝐫)}}{\displaystyle 𝑑𝐫^{}𝐚(𝐫^{})\psi (𝐫^{})\psi ^{}(𝐫^{})}`$
$`=`$ $`{\displaystyle \frac{1}{2i}}𝐚(𝐫)\psi (𝐫){\displaystyle \frac{1}{2i}}{\displaystyle \frac{\delta }{\delta \psi ^{}(𝐫)}}{\displaystyle 𝑑𝐫^{}\left\{𝐚(𝐫^{})\psi (𝐫^{})\psi ^{}(𝐫^{})\right\}}`$
$`+`$ $`{\displaystyle \frac{1}{2i}}{\displaystyle \frac{\delta }{\delta \psi ^{}(𝐫)}}{\displaystyle 𝑑𝐫^{}\psi ^{}(𝐫^{})\left\{𝐚(𝐫^{})\psi (𝐫^{})\right\}}`$
$`=`$ $`{\displaystyle \frac{1}{2i}}𝐚(𝐫)\psi (𝐫){\displaystyle \frac{1}{2i}}{\displaystyle \frac{\delta }{\delta \psi ^{}(𝐫)}}{\displaystyle 𝑑S𝐧𝐚(𝐫^{})\psi (𝐫^{})\psi ^{}(𝐫^{})}+{\displaystyle \frac{1}{2i}}\{𝐚(𝐫)\psi (𝐫)\}`$
$`=`$ $`{\displaystyle \frac{1}{2i}}𝐚(𝐫)\psi (𝐫)+{\displaystyle \frac{1}{2i}}\{𝐚(𝐫)\psi (𝐫)\}={\displaystyle \frac{1}{2i}}[𝐚(𝐫),]_+\psi (𝐫)`$
## Appendix B Matrix elements of the anti-commutator $`[𝐚(𝐫),]_+`$

The computation of the matrix element of the anti-commutator is presented below. The notation $`^a`$ means that the gradient acts on the Lagrange multiplier $`𝐚(𝐫)`$ only.

$`\varphi _n|[𝐚(𝐫),]_+|\varphi _m=\varphi _n|^a𝐚(𝐫)+2𝐚(𝐫)|\varphi _m`$
$`={\displaystyle 𝑑𝐫\varphi _n^{}(𝐫)\varphi _m(𝐫)𝐚(𝐫)}+2{\displaystyle 𝑑𝐫𝐚(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)}`$
$`={\displaystyle 𝑑𝐫\{𝐚(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)\}}{\displaystyle 𝑑𝐫𝐚(𝐫)\{\varphi _n^{}(𝐫)\varphi _m(𝐫)\}}+2{\displaystyle 𝑑𝐫𝐚(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)}`$
$`={\displaystyle 𝑑S𝐧𝐚(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)}{\displaystyle 𝑑𝐫𝐚(𝐫)\left(\varphi _m(𝐫)\varphi _n^{}(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)\right)}`$
$`={\displaystyle 𝑑𝐫𝐚(𝐫)\left(\varphi _m(𝐫)\varphi _n^{}(𝐫)\varphi _n^{}(𝐫)\varphi _m(𝐫)\right)}`$
Table 1. The energy of the active, i.e current carrying $`LxR`$, part of the system as a function of the steady current. $`\mathrm{\Delta }E_d=E_d(I)E_d(I=0)`$ is the energy increase to establish of the current $`I`$ carrying state.

| $`I`$ | $`E_d(I)`$ | $`\mathrm{\Delta }E_d(I)`$ |
| --- | --- | --- |
| 0.00 | 0.4977 | 0.0000 |
| 0.01 | 0.4994 | 0.0017 |
| 0.02 | 0.5043 | 0.0066 |
| 0.03 | 0.5124 | 0.0147 |
| 0.04 | 0.5233 | 0.0256 |
| 0.05 | 0.5368 | 0.0391 |
| 0.06 | 0.5526 | 0.0549 |
| 0.07 | 0.5706 | 0.0729 |
| 0.08 | 0.5905 | 0.0928 |
| 0.09 | 0.6122 | 0.1145 |
| 0.1 | 0.6356 | 0.1379 |

Figure 1. The Lagrange multiplier $`a(x)`$ is plotted as a function of the coordinate $`x`$. Calculations were performed for the step function current density distribution $`j=0.1\theta (x2,2x)`$.

Figure 2. The spatial dependences of real and imaginary components of the current carrying wave function. Calculations were performed for the step function current density distribution $`j=0.1\theta (x2,2x)`$. The dashed line on the upper panel is the ground state wave function of the system with zero current.

Figure 3. The spatial dependences of imaginary components of the wave functions are plotted for different values of the current. The current density is constrained to be the step function $`j=I\theta (x2,2x)`$. The solid line corresponds to $`I=0.01`$, dotted to $`I=0.05`$ and dashed to $`I=0.1`$.


text:
# On ghost condensation, mass generation and Abelian dominance in the Maximal Abelian Gauge

## I Introduction

As everybody knows, quarks are confined: nature as well as lattice simulations of nature are telling us that. Still, there is no rigorous proof of confinement. One proposal for the explanation of confinement is the idea of the dual superconductor: magnetic monopoles condense and induce a dual Meissner effect: color-electric flux between charges is squeezed and a string is created in between. The original work on this topic can be found in Nambu:1974zg ; Mandelstam:1974pi ; thooft . Abelian projection 'tHooft:1981ht  is a way to reveal the relevant degrees of freedom (the monopoles). In a lose way of speaking, at points were the projection is ill-defined, singularities invoke (Abelian) monopoles. Abelian dominance means that low energy QCD is dominated by Abelian degrees of freedom. Some early work on this is presented in Ezawa:bf . Numerical evidence can be found in e.g. Suzuki:1992gz ; Hioki:1991ai ; Suzuki:1989gp  and more recently Amemiya:zf .
Can this Abelian dominance be founded on more theoretical grounds? In the light of renormalization $`\stackrel{`}{\text{a}}`$ la Wilson, and assuming that the *off-diagonal* gluons (ghosts) attain a mass $`M`$ while the *diagonal* ones remain massless, an effective theory in terms of the massless diagonal fields could be achieved at low energy ($`M`$), thereby realizing a kind of Abelian dominance. In the context of low energy theories, we like to refer to the Appelquist-Carazzone decoupling theorem Appelquist:tg , which states that heavy particle modes decouple at low energy. Notice that this decoupling does not mean ”heavy terms” are simply removed by hand from the Lagrangian, their influence is still present through renormalization effects. As an illustration of this: a low energy, *Abelian* theory for Yang-Mills was derived in Kondo:2000sh , but the corresponding $`\beta `$-function was shown to be the same as the full Yang-Mills one.
The aforementioned pathway has been followed in a series of papers by Kondo et al Kondo:2000sh ; Kondo:2000zv ; Kondo:2000ey ; Kondo:1999tj ; Kondo:1999je ; Kondo:1997kn ; Kondo:1998sr ; Kondo:1997pc  and more recently the technique of the exact renormalization group has been employed by Freire Freire:2002ba ; Freire:2001nd  to construct effective low energy descriptions of Yang-Mills theory. The results have been used in order to construct a linearly rising potential between static quarks, a criterion for confinement. Their efforts were based on the dual superconductor picture, realized with MAG. Also the monopole condensation was discussed in their framework. An essential ingredient of their work is the mass scale of the off-diagonal fields. The monopole condensate is proportional to this mass squared Kondo:2000sh . The lattice reported a value of approximately 1.2 GeV for the off-diagonal gluon mass in MAG Yang-Mills Amemiya:zf . Next to these numerical results, analytical information is needed how this mass raises. A few papers have been written on this issue Kondo:2000ey ; Schaden:2001xu ; Schaden:2000fv ; Schaden:1999ew . All these authors came to the same conclusion: a dimension two ghost condensation gives an off-diagonal mass M. We already mentioned (but did not show explicitly) in a previous paper that we found the ghost condensation gives a *tachyonic* off-diagonal gluon mass Dudal:2002aj . In this paper, we will perform the calculations explicitly step by step. To make it self-contained, we will start from the beginning and in order to make comparison as transparent as possible, we will follow the (notational) conventions of Kondo:2000ey . For the sake of simplicity, we will restrict ourselves to the SU(2) case. We discuss the (in)completeness of presented work. We end with the path we intend to follow in the future to investigate dynamical mass generation in MAG.

## II Ghost condensation in the Maximal Abelian Gauge

Consider the Yang-Mills Lagrangian in four-dimensional Minkowski space time

$$=\frac{1}{4}_{\mu \nu }^A^{A\mu \nu }+_{GF+FP}$$
(1)
where $`_{GF+FP}`$ is the gauge fixing and Faddeev-Popov part.
We decompose the gauge field as

$`𝒜_\mu `$ $`=`$ $`𝒜_\mu ^AT^A=a_\mu T^3+A_\mu ^aT^a`$ (2)
$`_{\mu \nu }`$ $`=`$ $`_{\mu \nu }^AT^A=𝒟_\mu 𝒜_\nu 𝒟_\nu 𝒜_\mu =_\mu 𝒜_\nu _\nu 𝒜_\mu ig[𝒜_\mu ,𝒜_\nu ]`$ (3)
The $`T^A`$’s are the Hermitian generators of SU(2) and obey the commutation relations $`[T^A,T^B]=if^{ABC}T^C`$. $`T^3`$ is the diagonal generator. The capital index $`A`$ runs from 1 to 3. Small indices like $`a`$, $`b`$,… run from 1 to 2 and label the off-diagonal components. We will drop the index 3 later on.
As a gauge fixing procedure, we use MAG. Introducing the functional

$$[A]=(VT)^1d^4x\left(\frac{1}{2}A_\mu ^aA^{\mu a}\right)$$
(4)
with $`VT`$ the space time volume, MAG is defined as that gauge which minimizes $``$ under local gauge transformations. Since (4) is invariant under U(1) transformations w.r.t. the ”photon” $`a_\mu `$, MAG is only a partial gauge fixing. We do not fix the residual U(1) gauge freedom, since it plays no role for what we are discussing here.
To implement the gauge fixing in the Lagrangian (1), we use the so-called modified MAG. This gauge is slightly different from the ordinary MAG, it possesses for instance some more symmetry (see Kondo:1998sr  and references therein). Moreover, it generates the four-point ghost interaction, indispensable for the renormalizibility of MAG, as was proven in Min:bx .
Explicitly, we get

$$_{GF+FP}=i\delta _B\overline{\delta }_B\left(\frac{1}{2}A_\mu ^aA^{\mu a}\frac{\alpha }{2}iC^a\overline{C}^a\right)$$
(5)
where $`\alpha `$ is a gauge parameter, $`C`$ and $`\overline{C}`$ denote the (off-diagonal) ghosts and anti-ghosts, $`\delta _B`$ and $`\overline{\delta }_B`$ are the BRST and anti-BRST transformation respectively, defined by 111$`𝒞=(C^a,C^3)`$ with $`C^3`$ the diagonal ghost. Analogously for $`\overline{𝒞}`$.

$`\delta _B𝒜_\mu `$ $`=`$ $`𝒟_\mu 𝒞=_\mu 𝒞ig[𝒜_\mu ,𝒞]`$
$`\delta _B𝒞`$ $`=`$ $`{\displaystyle \frac{ig}{2}}[𝒞,𝒞]`$
$`\delta _B\overline{𝒞}`$ $`=`$ $`i`$
$`\delta _B`$ $`=`$ $`0`$ (6)
$`\overline{\delta }_B𝒜_\mu `$ $`=`$ $`𝒟_\mu \overline{𝒞}=_\mu \overline{𝒞}ig[𝒜_\mu ,\overline{𝒞}]`$
$`\overline{\delta }_B\overline{𝒞}`$ $`=`$ $`{\displaystyle \frac{ig}{2}}[\overline{𝒞},\overline{𝒞}]`$
$`\overline{\delta }_B𝒞`$ $`=`$ $`i\overline{}`$
$`\overline{\delta }_B\overline{}`$ $`=`$ $`0`$
$`+\overline{}`$ $`=`$ $`g[𝒞,\overline{𝒞}]`$ (7)
with the following properties

$`\delta _B^2`$ $`=`$ $`\overline{\delta }_B^2=\{\delta _B,\overline{\delta }_B\}=0`$
$`\delta _B\left(XY\right)`$ $`=`$ $`\delta _B\left(X\right)Y\pm X\delta _B\left(Y\right)`$
$`\overline{\delta }_B\left(XY\right)`$ $`=`$ $`\overline{\delta }_B\left(X\right)Y\pm X\overline{\delta }_B\left(Y\right)`$ (8)
where the upper sign is taken for bosonic $`X`$, and the lower sign for fermionic $`X`$.
Performing the BRST and anti-BRST transformations, yields

$`_{GF+FP}`$ $`=`$ $`B^aD_\mu ^{ab}A^{\mu b}+{\displaystyle \frac{\alpha }{2}}B^aB^a+i\overline{C}^aD_\mu ^{ac}D^{\mu cb}C^big^2ϵ^{ad}ϵ^{cb}\overline{C}^aC^bA^{\mu c}A_\mu ^d`$ (9)
$`+`$ $`i\overline{C}^agϵ^{ab}C^3D_\mu ^{bc}A^{\mu c}i\alpha gϵ^{ab}B^a\overline{C}^bC^3+{\displaystyle \frac{\alpha }{4}}g^2ϵ^{ab}ϵ^{cd}\overline{C}^a\overline{C}^bC^cC^d`$
where

$$D_\mu ^{ab}D_\mu ^{ab}[a]=_\mu \delta ^{ab}gϵ^{ab}a_\mu $$
(10)
is the covariant derivative w.r.t. the U(1) symmetry and

$`ϵ^{12}`$ $`=`$ $`ϵ^{21}=1`$ (11)
$`ϵ^{11}`$ $`=`$ $`ϵ^{22}=0`$ (12)
When we integrate the multipliers $``$ out, we finally obtain

$$_{GF+FP}=\frac{1}{2\alpha }\left(D_\mu ^{ab}A^{\mu b}\right)^2+i\overline{C}^aD_\mu ^{ac}D^{\mu cb}C^big^2ϵ^{ad}ϵ^{cb}\overline{C}^aC^bA^{\mu c}A_\mu ^d+\frac{\alpha }{4}g^2ϵ^{ab}ϵ^{cd}\overline{C}^a\overline{C}^bC^cC^d$$
(13)
Notice that the diagonal ghost $`C^3`$ has dropped out of (13).
For the (singular) choice $`\alpha =0`$, the 4-ghost interaction cancels from the Lagrangian. However, radiative corrections due to the other, non-vanishing 4-point interactions, reintroduce this term. We further assume that $`\alpha 0`$. Some more details concerning the properties for $`\alpha =0`$ can be found in Schaden:1999ew .
To discuss the ghost condensation mechanism, we ”Gaussianize” the 4-ghost interaction in the Lagrangian by means of the (U(1) invariant) auxiliary field $`\varphi `$

$$\frac{\alpha }{4}g^2ϵ^{ab}ϵ^{cd}\overline{C}^a\overline{C}^bC^cC^d\frac{1}{2\alpha g^2}\varphi ^2i\varphi ϵ^{ab}\overline{C}^aC^b$$
(14)
A useful identity to prove (14), reads

$$ϵ^{ab}ϵ^{cd}\overline{C}^a\overline{C}^bC^cC^d=2\left(iϵ^{ab}\overline{C}^aC^b\right)^2$$
(15)
The part of the Lagrangian which concerns us for the moment is

$$\widehat{}=i\overline{C}^a_\mu ^\mu C^a\frac{1}{2\alpha g^2}\varphi ^2i\varphi ϵ^{ab}\overline{C}^aC^b$$
(16)
Assuming constant $`\varphi `$, we use the Coleman-Weinberg construction Coleman:jx  of the effective potential $`V(\varphi )`$. This means we are summing all 1-loop (off-diagonal) ghost bubbles with any number of $`\varphi `$-insertions. This yields

$$(VT)V(\varphi )=d^4x\frac{\varphi ^2}{2\alpha g^2}+i\mathrm{ln}det\left(_\mu ^\mu \delta ^{ab}\varphi ϵ^{ab}\right)$$
(17)
or

$$V(\varphi )=\frac{\varphi ^2}{2\alpha g^2}\frac{1}{i}\frac{d^4k}{\left(2\pi \right)^4}\mathrm{ln}\left(k^4+\varphi ^2\right)$$
(18)
Employing the Wick rotation $`k_0ik_0`$ 222If one would like to avoid Wick rotations, one could start immediately from the Euclidean version of Yang-Mills., and performing the integration in dimensional regularization within the $`\overline{MS}`$ scheme, we arrive at

$$V(\varphi )=\frac{\varphi ^2}{2\alpha \overline{g}^2}+\frac{\varphi ^2}{32\pi ^2}\left(\mathrm{ln}\frac{\varphi ^2}{\overline{\mu }^4}3\right)$$
(19)
This potential possesses a local maximum at $`\varphi =0`$ (the usual vacuum), but has global minima at

$$\varphi =\pm v=\pm \overline{\mu }^2e^{1\frac{8\pi ^2}{\alpha \overline{g}^2(\overline{\mu }^2)}}$$
(20)
We take $`\alpha >0`$ since $`v`$ diverges for $`\overline{g}^20`$ if $`\alpha <0`$.
Up to now, we find complete agreement with Kondo:2000ey . We proceed by calculating the ghost propagator in the non-zero vacuum ($`V(v)<0`$). Substituting $`\varphi =v`$ in (16), it is straightforward to determine the Feynman propagator

$$C^a(x)\overline{C}^b(y)=\frac{d^4k}{\left(2\pi \right)^4}\frac{k^2\delta ^{ab}+vϵ^{ab}}{k^4+v^2}e^{ik(xy)}$$
(21)
With the above propagator, we are ready to determine the 1-loop off-diagonal gauge boson polarization. Now, there exists a non-trivial contribution coming from the ghost bubble, originating in the interaction term $`ig^2ϵ^{ad}ϵ^{cb}\overline{C}^aC^bA^{\mu c}A_\mu ^d`$, resulting in a mass $`M`$ for the off-diagonal gluons. Again Wick rotating $`k_0ik_0`$ to get an integral over Euclidean space time, one easily obtains

$$M^2=g^2\frac{d^4k}{\left(2\pi \right)^4}\frac{2k^2}{k^4+v^2}$$
(22)
There is one remaining step, we still have to calculate the integral of (22). Using dimensional regularization, we find the finite result

$$M^2=\frac{g^2v}{16\pi }<0$$
(23)
where we have used that $`v>0`$. Here we find a *different* result in comparison with the other references Kondo:2000ey ; Schaden:2001xu ; Schaden:2000fv ; Schaden:1999ew . To be more precise, we find the opposite sign. This sign difference is not meaningless, since the negative sign we find means that the off-diagonal fields have a *tachyonic* mass.
Hence, we state that a ghost condensation $`\stackrel{`}{\text{a}}`$ la $`ϵ^{ab}\overline{C}^aC^b`$ is *not* the mechanism behind the off-diagonal mass generation in MAG, and consequently does not give evidence for Abelian dominance.
Another important point is what happens with the diagonal gluon. Consider the term $`i\overline{C}^aD_\mu ^{ac}D^{\mu cb}C^b`$ of (13), it contains a part proportional to $`i\overline{C}^aC^aa_\mu a^\mu `$. Doing the same as for the off-diagonal gluons, the diagonal gluon $`a_\mu `$ seems to get a (real) mass too, which is of the same order as the off-diagonal one (up to the sign). However, there are other 1-loop contributions coming from the terms proportional to $`ϵ^{ab}\left(_\mu \overline{C}^a\right)C^ba^\mu `$ and $`ϵ^{ab}\overline{C}^a\left(_\mu C^b\right)a^\mu `$. These contributions cancel the one coming from the term proportional to $`i\overline{C}^aC^aa_\mu a^\mu `$. Consequently, the ”photon” $`a_\mu `$ remains massless, as could be expected by the residual U(1) invariance.
Another point of concern is the renormalizibility of the ”Gaussianized” Lagrangian. A completely analogous approach can be done in case of the 2-dimensional Gross-Neveu model Gross:jv , where the 4-fermion interaction can also be made Gaussian by the introduction of an auxiliary field $`\sigma `$. This works well at 1-loop order, but from 2 loops on, *ad hoc* counterterms have to be added in order to end up with finite results Luperini:1991sv . A successful formalism to deal with local composite operators in case of the Gross-Neveu model was developed in Verschelde:jx . A similar approach should be used to investigate the ghost condensates.
One could wonder what the mechanism behind the mass generation might be, since the previous paragraphs showed that we didn’t find a dynamically generated real mass for the (off-diagonal) particles. In order to find an answer to this question, we first give a very short overview of recent results in the Landau gauge, giving us a hint in which direction we should look for the mass generation.

## III  Gluon condensation via $`𝒜^2`$ in the Landau gauge and its nephew $`A^2`$ in the Maximal Abelian Gauge

A well known condensate in QCD (or Yang-Mills) is the dimension four gluon condensate $`_{\mu \nu }^𝒜^{\mu \nu 𝒜}`$. This is the lowest dimensional gluonic condensate that can exist, since no local, gauge-invariant condensates with dimension lower than 4 exist. However, recently interest arised concerning a dimension 2 gluon condensate in Yang-Mills theory in the Landau gauge. One way it came to attention was the conclusion that there exists a non-negligible discrepancy between the lattice strong coupling constant $`\alpha _s`$ (determined via the 3-point gluon interaction) and the perturbative one, into a relatively high energy region where this wouldn’t be expected (up to 10 GeV). Also the propagator showed a similar discrepancy. The $`\frac{1}{p^4}`$ power correction due to $`^2`$ is far to small to explain this. It was shown that a $`\frac{1}{p^2}`$ power correction could solve the discrepancy. More precisely, the Operator Product Expansion (OPE) used in combination with the $`𝒜^2`$ condensate was able to fit both predictions Boucaud:2000 ; Boucaud:2001st . An important question that naturally arises, sounds: has $`𝒜^2`$ any physical meaning, or is it merely a gauge artefact? The point is that $`𝒜^2`$ equals $`(VT)^1\mathrm{min}_Ud^4x𝒜_\mu ^U𝒜^{\mu U}`$ in the Landau gauge, and this latter operator is, although non-local, gauge-invariant. Hence, $`𝒜^2`$ can be given some physical sense in the Landau gauge. Moreover, Gubarev:2000eu  discussed the relevance of $`𝒜^2`$ in connection with topological structure (monopoles) of compact QED. The physical relevance of the Landau gauge, in the framework of geometrical monopoles, is explained in Gubarev:2000nz . The authors of that paper also stress that the values found with an OPE calculation, only describe the soft (infrared) content of $`𝒜^2`$, while they argue that also hard (short range) contributions, unaccessible for OPE, may occur. In this context, we cite Verschelde:2001ia , where a formalism was constructed for the calculation of the vacuum expectation value of (local) *composite* operators. Since this is based on the effective action, it should in principle, give the ”full” value of $`𝒜^2`$, i.e. soft and hard part. For example, one could assume an instanton background as a possible source of long range contributions. In fact, there is some preliminary evidence that instantons can explain the OPE values Boucaud:2002nc .
The conclusion that one can draw from all this is that the dimension 2 condensate $`𝒜^2`$ may have some physical relevance in the Landau gauge.
Let us go back to MAG 333Here, with MAG we mean the gauge minimizing the functional (4), and not the modified MAG. now. In this particular gauge, $`(VT)^1\mathrm{min}_Ud^4x𝒜_\mu ^U𝒜^{\mu U}`$ no longer reduces to a local operator. It would be interesting to repeat e.g. the OPE calculations of Boucaud:2001st  for the coupling constant and propagators in MAG, but which dimension 2 condensate(s) could take over the role of $`𝒜^2`$ in the Landau gauge? To solve this, we draw attention to the striking similarity existing between the Landau gauge and MAG. The former one can be seen as that gauge minimizing $`(VT)^1d^4x\left(𝒜_\mu ^A𝒜^{\mu A}\right)^U`$, while the latter one minimizes $`(VT)^1d^4x\left(A_\mu ^aA^{\mu a}\right)^U`$. This operator reduces to the local one $`A^2`$ in MAG and can be seen as the MAG version of $`𝒜^2`$. Due to the more complex nature of the (renormalizable) modified MAG, other dimension 2 condensates exist (the ghost condensates). Notice that all these condensates are U(1) invariants, hence the U(1) symmetry will be preserved.
The physics we see behind all these condensates is that they might have a common, deeper reason for existence. In this context, we quote Fukuda:1977mz ; Fukuda:1977zp , where it was shown that the zero vacuum is instable (tachyonic) and a vacuum with lower energy is achieved through gluon pairing, and an accompanying gluon mass. The vacuum energy itself is a physical object. After choosing a certain gauge, the different types of dimension 2 condensates are just an expression of the fact that $`E=0`$ is a wrong vacuum state. In this sense, all these different condensates in different gauges are equivalent in a way, since they lower the vacuum energy to a stable $`E<0`$ vacuum. In Gubarev:2000nz , discussion can be found on the appearance of the soft part of $`𝒜^2`$ in gauge-variant quantities like in an OPE improvement of the gluon propagator, while the hard part enters physical quantities. It is imaginable that the mechanism behind this hard part (see Gubarev:2000nz  for more details) is the same in different gauges, but reveals its importance with different condensates, depending on the specific gauge. This might justify its possible appearance in gauge-invariant quantities.

## IV  Further discussion on the ghost condensation and mass generation in the modified MAG

In Dudal:2002ye , it was shown that it is possible to fix the residual abelian gauge freedom of MAG in such a way that the ghost condensate $`ϵ^{ab}\overline{C}^aC^b`$ does not give rise to any mass term. This abelian gauge fixing (needed for a complete quantization of the theory) was based on the requirement that the fully gauge fixed Lagrangian has a SL(2,R) and anti-BRST invariance. A restricted 444By restricted, we mean that the symmetry only acts non-trivially on the off-diagonal fields.version of this SL(2,R) symmetry was originally observed in SU(2) MAG in Schaden:1999ew , and later generalized to SU(N) MAG Sawayanagi:zw . In Dudal:2002ye , the symmetry was defined on all the fields (diagonal and off-diagonal). In fact, that SL(2,R) symmetry together with the (anti-) BRST symmetry form a larger algebra, the Nakanishi-Ojima (NO) algebra. This NO algebra is known to generate a symmetry of the Landau gauge and a certain class of generalized covariant gauges, more precisely the Curci-Ferrari gauges, given by the gauge fixing Lagrangian

$$_{GF+FP}=i\delta _B\overline{\delta }_B\left(\frac{1}{2}𝒜_\mu ^A𝒜^{\mu A}\frac{\alpha }{2}i𝒞^A\overline{𝒞}^A\right)$$
(24)
The Landau gauge corresponds to the gauge parameter choice $`\alpha =0`$. For more details, see no ; ni ; oji ; dj ; ds .
Yang-Mills theory with the gauge fixing (24) possesses a generalization to a massive SU(N) gauge model, the so-called Curci-Ferrari model cf . Although this model is non-unitary, it is known to be (anti-)BRST invariant and renormalizible, whereby the mass term is of the form

$$_{mass}=M^2\left(\frac{1}{2}𝒜_\mu ^A𝒜^{\mu A}i\alpha 𝒞^A\overline{𝒞}^A\right)$$
(25)
Keeping this in mind and recalling that in Verschelde:2001ia , a dynamically generated mass was found in case of the Landau gauge by coupling a source $`J`$ to the operator $`A^2`$, it becomes clear that in case of the Curci-Ferrari gauge, the same technique could be employed by coupling a source $`J`$ to the composite operator

$$_{source}=J\left(\frac{1}{2}𝒜_\mu ^A𝒜^{\mu A}i\alpha 𝒞^A\overline{𝒞}^A\right)$$
(26)
Returning to the case of MAG and comparing the gauge fixing Lagrangians (5) and (24), the equivalent of (26) reads

$$_{source}=J\left(\frac{1}{2}A_\mu ^aA^{\mu a}i\alpha C^a\overline{C}^a\right)$$
(27)
This idea to arrive at a dynamically generated mass in case of the Curci-Ferrari and Maximal Abelian gauge was already proposed in Kondo:2001nq ; Kondo:2001tm . There, it was explicitly shown that the operator coupled to the source $`J`$ in the expressions (26) or (27), is on-shell BRST invariant.
We reserve the actual discussion of the aforementioned framework to get a dynamical mass for future publications, since it is quite involved and a clean treatment of it needs a combination of the local composite operator formalism Verschelde:2001ia  and the algebraic renormalization technique book ; bbh .
Before turning to conclusions, we want to draw attention to the following. We decomposed the 4-ghost interaction with a real auxiliary field $`\varphi `$ whereby $`\varphi ϵ^{ab}\overline{C}^aC^b`$. Let’s make a small comparison with ordinary superconductivity. Usually, there is talked about BCS pairing, i.e. particle-particle and hole-hole pairing. The analogy of this in the ghost condensation case would be ghost-ghost pairing and antighost-antighost pairing. This can be achieved by an alternative decomposition of the 4-ghost interaction via a pair of auxiliary fields $`\sigma `$ and $`\overline{\sigma }`$ such that $`\sigma ϵ^{ab}C^aC^b`$ and $`\overline{\sigma }ϵ^{ab}\overline{C}^a\overline{C}^b`$. This kind of pairing 555Our conclusion about the tachyonic mass is unaltered by this alternative decomposition of the 4-ghost interaction. was considered in Lemes:2002ey . A less known effect is the particle-hole pairing, the so-called Overhauser pairing bo . This corresponds to the kind of condensation we and the papers Kondo:2000ey ; Schaden:1999ew ; Schaden:2000fv ; Schaden:2001xu  considered. From the viewpoint of the SL(2,R) symmetry, the existence of different channels where the ghost condensation can take place should not be suprising. The different composite ghost operators are mutually changed into each other under the action of the symmetry. Here and in the other papers the *choice* was made to work with the Overhauser channel, but a complete treatment would need an analysis of all channels at once, and with the local composite operator technique. This analysis of the BCS versus Overhauser effect is nicely intertwined with the existence of the NO algebra and its (partial) breakdown, and it is very much alike for the MAG, Landau bezig  and Curci-Ferrari gauge, just as in case of the mass generation mechanism. As an indication, it has been found recently that, although no 4-ghost interaction is present in the Landau gauge, the condensation $`\stackrel{`}{\text{a}}`$ la $`f^{ABC}\overline{𝒞}^A𝒞^B`$ etc. also occurs Lemes:2002rc .

## V  Conclusion

We considered Yang-Mills theory in the Maximal Abelian Gauge. With this non-linear gauge choice, a 4-ghost interaction enters the Lagrangian. Such an interaction could allow a non-zero vacuum expectation value for (off-diagonal) dimension 2 ghost condensates. Consequently, it was expected that a mass generating mechanism for the off-diagonal gluons *and* the diagonal gluons due to 4-point interaction terms of the form gluon-gluon-ghost-anti-ghost was found.
We explained why this particular type of ghost condensation is not sufficient to construct a (off-diagonal) dynamical mass in SU(2) Yang-Mills theory in the Maximal Abelian Gauge, an indicator for Abelian dominance. We have restricted ourselves to the SU(2) case, but a similar conclusion will exist for general SU($`N`$). Explicit calculations showed that we ended up with a tachyonic off-diagonal mass $`M`$ ($`M^2<0`$). This result indicate something is missing. A comparison with Yang-Mills theory in the Landau gauge and the role played by the mass dimension 2 gluon condensate $`𝒜^2`$, shed some light on the route that should be followed.
We revealed certain shortcomings of the present available studies on the ghost condensation (renormalizibility, existence of more than one condensation channel).
The actual study of the mass generation and the ghost condensation with its symmetry breaking pattern will be discussed elsewhere. We will follow the local composite operator formalism of Verschelde:2001ia , where a source is coupled to each operator and the effective action can be treated consistently. This effective potential formalism allows a clean treatment of the role played by the dimension 2 operators. We remark that with essentially perturbative techniques one can obtain at least qualitatively trustworthy results 666In case of the Gross-Neveu model, very accurate results were obtained Verschelde:jx . In case of the $`A^2`$ condensate in the Landau gauge, the relevant coupling constant was quite small, making the expansion acceptable Verschelde:2001ia . on the stability of the condensates and their relevance for e.g. mass generation and symmetry breakdown, without making it directly necessary to go to (or extrapolating to) strong coupling.
We conclude by mentioning that the dimension 2 condensates and the accompanying mass generation in Yang-Mills are not only of theoretical importance (the role of $`𝒜^2`$ for OPE corrections Boucaud:2000 ; Boucaud:2001st , monopoles Gubarev:2000eu ; Gubarev:2000nz , *short range linear correction* to the Coulomb-like potential Boucaud:2000 , low energy effective theories Kondo:2002xn ,…) but also have their importance for automated Feynmandiagram calculations Gracey:2001 ; Browne:2002wd ; Gracey:2002  where a gluon mass serves as a infrared regulator. If this mass is generated in massless Yang-Mills, it does not have to be implemented by hand.

## VI Acknowledgments

Private communication with M. Schaden, S. P. Sorella and J. A. Gracey is greatly acknowledged.