\section{Introduction}

X-ray surveys have
consistently indicated the existence of a moderately bright X-ray
source near the galactic longitude 41.97$^{\circ}$ and
latitude $-0.80^{\circ}$.  This source, which we refer to as X 1908+075,  was cataloged as 4U 1909+07 by Uhuru, the first earth-orbiting
mission dedicated to celestial X-ray astronomy
(\markcite{forman78}Forman et al. 1978). Subsequent determinations were cataloged as 1M 1912+077 with MIT/OSO 7 (\markcite{markert79}Markert et al. 1979), 3A 1907+074 with Ariel 5 (\markcite{warwick81}Warwick et al. 1981; \markcite{bell84}Bell Burnell \& Chiappetti 1984), 1H 1907+074 with HEAO A-1/LASS (\markcite{wood84}Wood et al. 1984), GPS
1908+075 with EXOSAT/ME (\markcite{warwick88}Warwick et al. 1988).  There are large uncertainties in the position determinations of these survey instruments, as shown in Fig. \ref{pos}.
The HEAO A-2 experiment on HEAO-1 scanned this region and detected X-ray emission of $\sim 13$  mCrab in the 1.5--20 keV band (HEASARC archival data).  The source was also detected with the modulation collimator experiment (HEAO A-3) on HEAO-1, with an 8$\sigma$ detection in each of the two collimators for the energy range of 3--13 keV (HEASARC archival data).  The HEAO A-3 instrument  yields a multiplicity of relatively precise error regions which appear as a grid of diamonds in Fig \ref{pos}.

The
Einstein imaging proportional counter (IPC) imaged
a portion of this region and detected a  source, 1E 1908.4+0730, the brightest source in the center of the field, at a significance level of $\sim  7\sigma $
 (HEASARC archival data). The IPC flux is about 0.3 mCrab (0.5--3.5 keV), possibly highly attenuated by interstellar absorption (see table 1).  The position (J2000) is $ \alpha =19^h 10^m 46^s$ and $\delta =+07^\circ 36'
07''$ with a positional accuracy of $50''$.    This position overlaps one of the HEAO A-3 ``diamonds''.   There are a few other 5-$\sigma$ detections at the edge of the field of view ($1^\circ$ square) of the IPC, and these sources are located well beyond the celestial map provided in Fig. \ref{pos}. Furthermore,  none of these sources have positions that are consistent with the HEAO A-3 diamonds.  It is thus highly likely that the Einstein/IPC detection represents the
same source as the X-ray survey detections.  No optical
counterpart has been identified within the error circle of 1E 1908.4+0730.

\setcounter{figure}{0}
\refstepcounter{figure}
\PSbox{fig1.eps hoffset=-15 voffset=-5 hscale=45 vscale=45}{3.0in}{3.0in}{\\\\\small Fig.1  Positions of X 1908+075 and other X-ray sources detected previously, with reported error boxes superposed on an optical image taken from  the Digitized Sky Survey.  The degenerated positions determined with HEAO A-3 can be seen as a grid of diamonds.  The position of X 1908+075 used in the ASM analysis was taken to be that of Einstein/IPC source 1E 1908.4+0730, represented as a small circle.
\label{pos}}

\vspace{0.5cm}

The X-ray spectra and flux available from the hard X-ray survey detections are
 roughly consistent with each other as summarized in table 1. The fluxes are about 5 mCrab in the 2--6 keV band and 7--14 mCrab in the 2--20 keV band. The spectrum of EXOSAT/ME observation of
 4U 1909+07 in the $1$--$15$  keV band (at the UHURU position) indicates excessive absorption
 at photon energy below $6$ keV compared with that of the Crab Nebula (HEASARC archival data). The
 observations from the Ariel 5 spectrometer show that 3A 1907+074 has a hard spectrum
 which can be fitted with a power law of photon index $1.6 \pm 0.15$,
 hydrogen column density of $3^{+1.5} _{-1.2} \times 10^{22}$
 cm$^{-2}$, and an iron line emission with equivalent width $660 \pm
 220$ eV (\markcite{bell84}Bell Burnell \& Chiappetti 1984). None of these pointed observations were long enough to investigate
 the long-term variability of this source.

In this paper, we use the data from the All Sky Monitor (ASM) (\markcite{levine96}Levine et al. 1996) on board {\it RXTE} (\markcite{bradt93}Bradt, Rothschild, \& Swank 1993) to confirm the position and to study the nature of X 1908+075. We show that the coded mask patterns for the ASM,  processed for the
 Einstein/IPC position, yield a persistent $\sim 8$ mCrab X-ray
 source with a coherent $4.4$-day periodicity for over 3
 years. We then discuss the implications of
 this discovery.

\section{Observation}
The {\it RXTE}/ASM has been monitoring the sky
routinely since $1996$ February.  It consists of three Scanning Shadow Cameras, each
consisting of a coded mask and a position-sensitive proportional
counter. The ASM is well-suited for our purpose since each of its cameras  has a 1.5--12 keV response, an intrinsic angular resolution of a few arc minutes, and a large field of view of $6^\circ \times 90^\circ$(FWHM). In addition,   the ASM  provides observations with
relatively good sensitivity ($\le 5$--$10 $ mCrab in a day), frequent
data sampling ($\sim 5$--30 times a day) and long baseline (thus far 3
years).

For a 90-s exposure in each ASM camera, a set of source
intensities for the known sources in the FOV is derived from a $\chi ^2$ fit for the overlapping
shadow patterns of the X-ray sources and the X-ray background.  The intensity is computed in each of the three energy
bands ($1.5$--$3$, $3$--$5$ and $5$--$12$ keV) and the results are
normalized in units of the count rate for a source at the center of
the field of view of camera \# 1. In these units, the
$1.5$--$12$ keV flux for the Crab Nebula is $75.5$ ASM ct/s. The
estimated errors of the source intensities include the uncertainties
due to counting statistics and a systematic error of $1.9\%$
determined empirically from the light curves of the Crab Nebula. A
detailed description of the ASM and the light curves can be found in
Levine et al. (1996) and Levine (1998).

\section{Analysis and Results}

The X-ray light curve from the ASM observations of X 1908+075 (1996 March to 1999 January) is shown in Fig. \ref{ltc_all}.  The intensity points are 5.0-d average of the intensities derived from the individual 90-s dwells for the assumed position of the Einstein/IPC source.  It is apparent that X 1908+075 has generally
persisted at an intensity level of about 0.6 ct/s, or 8.0 mCrab in the 1.5-12 keV
band for about three  years.  The spatial  response function of the ASM implies that these photons arise within 12 arcmin of the Einstein/IPC source position, thus confirming its association with X 1908+075. The average X-ray intensity of X 1908+075 is well below the detection threshold (30 mCrab) for a single camera exposure  of the ASM. The effort to determine an independent source position using ASM data alone, e.g., with the superposition of cross-correlation maps,  would be limited by the appearance of systematic  noise feature with amplitudes that are a substantial fraction of the average source intensity.

The average values of the hardness ratios of X 1908+075 are HR1=$1.0 \pm 0.2$ , HR2=$3.6 \pm 0.3 $ respectively, where HR1 is defined as the ratio of the ASM count rates in the $3$--$5$ keV band to that of the
$1.5$--$3$ keV band, and HR2 is the ratio of the count rates in the
$5$--$12$ keV band relative to the $3$--$5$ keV band.
 X 1908+075 has a much  harder X-ray spectrum than either the Crab Nebula (HR1=0.9, HR2=1.1) or Cyg X-1 in the low-hard state (HR1=1.0, HR2=1.5).

\refstepcounter{figure}
\PSbox{fig2.eps hoffset=-10 voffset=0 hscale=52
vscale=52}{3.2in}{3.2in}{\\\small Fig. 2 $\it RXTE$/ASM light curves of X 1908+075 for the
period $1996$ February to $1998$ November.  The data are binned in $5.0$-day bins.  Only $5\sigma$ detected data points are plotted. MJD $51000$ (MJD=JD-2400000.5) corresponds to
$1998$ July $6$. There is a persistent X-ray emission at $\sim 0.6$ ct/s (8 mCrab) for  3 years. \label{ltc_all}}
\vspace{0.5cm}

Periodicities of X 1908+075 were sought by means of the Lomb-Scargle
periodograms of both the light curves and the derived hardness ratios.
The Lomb-Scargle periodogram (\markcite{lomb76}Lomb 1976; \markcite{scargle82}Scargle 1982; \markcite{press92}Press et al. 1992) was used for the calculation instead of the Fast
Fourier Transform (FFT) since the ASM data points are unevenly spaced
in time.  With this periodogram, a maximum in the power occurs at the
frequency which gives the least squares fit to a sinusoidal function. The expected average power for random noise is $1.0$. We
oversampled the power spectrum so that the frequencies are a factor of
4 more closely spaced than $1/T$, where $T$ is the total duration of
the data used.  The goal is to ensure the detection of a peak for a
signal that is of marginal statistical significance and to best
locate the peak.  The result is shown in Fig. \ref{power}.

There is a distinct narrow peak in the power spectra at a frequency
that is corresponding to $4.400$ d.  The width of the power peak is
0.009-d, consistent with that of a pure sinusoidal wave. The
uncertainty of this period is estimated to be $0.001$ d, assuming a
single signal with Gaussian noise. With 3 years of source coverage,
the uneven data spacing would not degrade the uncertainty to any
noticeable degree (\markcite{horne86}Horne \& Balinuas 1986).  If we assume the noise follows an
exponential distribution as expected for random noise, the false-alarm
probability of this detection is less than $ 10^{-7}$. In addition,
there is a low-frequency peak corresponding to about $180$-d
period. Except for these two peaks the overall power is flat with an
average of about $1.06$.

We have also conducted another independent periodicity search using
the epoch folding technique (\markcite{leahy83}Leahy et al. 1983; \markcite{davis90}Davies 1990).  The
data were folded modulo a trial period and then grouped by phases.  The resulting deviation from uniform random noise can
then be detected using a $\chi ^2 $ test.  For periods between 0.01--7
days, a distinct peak indicating the maximum deviation was found
around $4.400$ d with $< 10^{-6}$ false alarm probability level. This
is consistent with the results of the Lomb-Scargle periodogram.

\refstepcounter{figure}
\PSbox{fig3.eps  hoffset=-15 voffset=-30 hscale=55 vscale=60}{3.0in}{3.0in}{\\\small Fig. 3 The Lomb-Scargle periodogram of the light curves of X 1908+075. The dotted line indicates the $99.99999\%$ confidence level assuming an exponential distribution of the noise powers. There is a distinct peak at $4.4$-d period. The 180-d period is probably caused by scattered solar X-ray contamination. \label{power}}

\vspace{0.5cm}

The data were then folded modulo the period of the best-fit value of
$4.39986$ days from the Lomb-Scargle periodogram to study the
intensity variations with phase (Fig. \ref{fold}).

\refstepcounter{figure}
\PSbox{fig4.eps hoffset=-18 voffset=-15 hscale=50
vscale=50}{3.0in}{3.0in} {\\\\\small Fig. 4 Folded light curves of X 1908+075.  Observational data points are represented as histograms. Phase zero is defined as the phase of the minimum X-ray flux. The error bars represent one
standard deviation. The
solid lines are the best-fit sinusoidal waves. The most distinctive feature is the nearly sinusoidal modulation with large amplitude. \label{fold}}

\vspace{0.5cm}

 The most
distinctive feature in the folded light curves is the nearly  sinusoidal
modulation of large amplitude.  We fitted a sinusoidal wave plus a constant to the
folded light curves. The best fit yields a $26 \pm 3 \%$ amplitude of
the sinusoidal wave in the $1.5$--$12$ keV band, $51 \pm 14\%$ in the
$1.5$--$3$ keV band, $49 \pm 8\%$ in the $3$--$5$ keV band, and $22
\pm 3\%$ in the $5$--$12$ keV band. The best-fit epoch for the minimum
flux is:

\begin{center}
	$T_{0}=2450440.419 \pm 0.068 +n\times (4.400 \pm 0.001$ ) (JD)
\end{center}

The persistence and coherence of this periodic variation were tested
by subdividing the ASM light curves into four
consecutive $250$-d intervals and then searching for periodicities.  The same analyses described above were
applied to each of the $250$-d light curves. In the calculated power
density spectra, the signals are still visible as narrow peaks near
$P=4.40$ d in the power density spectra even though they are no longer
independently significant. The best-fit amplitudes and phases of the
sinusoidal waves are roughly consistent with each other (see Table
\ref{tab2} for a summary).  We thus conclude that the modulation is
persistent and coherent over the 3-year time interval.

\section{Discussion}

The X-ray emission detected by the ASM has been found to arise from the immediate region of the Einstein/IPC source, to exhibit a hard spectrum and intensity similar to the historical survey detections,  and to be persistent over $\sim 3$ years.  This together with the Einstein IPC and  HEAO A-3 positional correspondence strongly suggest that the position of 1E 1908.4+0730 is the
accurate location of the sources seen in the historical X-ray surveys.

It is natural to link the observed $4.4$-day period of X 1908+075 to
the orbital period of a binary. The periods of X-rays from known
sources are believed to be associated with one of the following
(\markcite{white95}White, Nagase, \& Parmar 1995): (1) rotational periods of neutron stars in X-ray
pulsars, ranging between milliseconds to hundreds of seconds; (2)
orbital periods of binary systems, ranging from tens of minutes to
tens of days, and possibly hundreds of days; and (3) ``superorbital''
orbital periods from $\sim 30$ days to hundreds of days, which in some cases are believed to be the precession period of the accretion
disk (e.g., \markcite{levine82}Levine \& Jernigan 1982).  It is clear that the coherence and the timescale of the $4.4$-day X-ray oscillation in X 1908+075 can be  best interpreted as the  orbital period of a binary system.

We can then picture X 1908+075 as a binary system consisting of a
compact object and a mass-donor companion star. An X-ray binary with a
low mass companion (LMXB) is unlikely the scenario for X 1908+075,
since the known orbital periods of most LMXBs (White et al. 1995) are on
the order of hours, much shorter than the detected 4.4 days. There are
only a few known LMXBs having orbital periods longer than a day, among
them only Cir X-1 ($P \sim 16.6$ d) and Her X-1 ($P \sim 1.7$ d, a
possible LMXB) exhibit orbital modulations in X-rays.  The modulation
in X-rays of Cir X-1 manifests itself as periodic flares and dips
likely due to a highly eccentric orbit, while Her X-1 is an eclipsing
system.  Such dramatic modulations are not evident in the orbital
light curves of X 1908+075. On the other hand, most HMXBs have orbital
periods on the orders of days.  The smooth and nearly sinusoidal
modulations in X 1908+075 have a lot in common with a supergiant
systems such as Cyg X-1.  We therefore suggest that X 1908+075 is a
high mass X-ray binary (HMXB).

We further suggest that the companion star is a supergiant rather than
a Be star for the following reasons.  (1) The $4.4$-d orbital period
fits in the range favored more by supergiant systems than by Be-star
systems, as the known orbital periods of the latter are much
longer than $4.4$ days ($\ge 16.7$-d) (White et al. 1995). (2) The known orbital light
curves of Be-star systems are all outburst-like rather than sinusoidal
(e.g., EXO 2030+375), probably due to the enhanced accretion caused by
the the compact object's periastron passage around the Be star.  (3) The orbital modulation of a supergiant system are known to be persistent while this is generally not true in
Be-star systems (e.g., 4U 0115+63 in ASM light curves during 1996 March
to 1998 November).

The smooth and roughly sinusoidal orbital modulation in the X-rays of
X 1908+075 suggests that the system could be embedded in a strong
stellar wind from the companion star (e.g., a supergiant O or B star),
which absorbs and scatters the X-rays. The X-ray modulation is then
caused by changes in the optical depth along the line of sight to the
X-ray source as a function of orbital phase.  For X 1908+075,
absorption of X-rays probably dominates over scattering, as the
amplitude of the modulation in the $1.5$--$5$ keV band is larger than in the $5$--$12$ keV band.  To explain a modulation with an amplitude of  nearly $\sim 22\%$ in the 5--12 keV band, the difference in
the hydrogen column densities between phases 0 and 0.5 is expected
to be more than  $10^{23}$ cm$^{-2}$ for an absorption and scattering process.
The mechanism of the orbital modulation of X 1908+075 could be very
similar to that of Cyg X-1, wherein   hard state orbital modulation in the ASM light
curves can be well explained with
absorption and scattering in a stellar wind (\markcite{wen99}Wen et al. 1999).  Another
example of such a system could be the newly discovered X-ray pulsar
XTE~ J1855-026 (\markcite{corbet99}Corbet et al. 1999).  In this system, the orbital period
is proposed to be 6.1-d as detected by the ASM.  The modulation is
roughly sinusoidal with an amplitude about $22\%$ in the 1.5--12 keV band.  It has
also been proposed as a supergiant wind accretion system.

The hydrogen column density $n_{H}$ can be used to estimate the
optical extinction in the V band using the empirical relation between
the $n_{H}$ and the color excess $E_{B-V}$ (\markcite{savage79}Savage \& Mathis 1979).  It is
found that $A_{V}=15$ for $n_{H} =3\times 10^{22}$ cm$^{-2}$, as
indicated from the Ariel 5 observation of 3A 1907+074. The optical
objects within the Einstein/IPC error circle have magnitudes $V \ge
19$ in the Palomar sky survey.  To estimate the distance to X 1908+075,
we first investigated the possible range of the absolute V-band
magnitudes $M_V$ (from \markcite{lang92}Lang 1992) known for supergiant stars in
binaries with known orbital periods.  We found that $M_V \le -5.2$.
Together with $A_{V}=15$ and $V \ge 19$, we obtained $0.7$ kpc as the
lower-limit for the distance. The X-ray luminosity, estimated from the
ASM count rate (1.5--12 keV), is then $L_{1.5-100 keV} \ge 5 \times
10^{34}$ erg s$^{-1}$ if we assume the Ariel 5 spectrum. Given the
$A_V$ value, a more likely scenario is a distance $\sim$ 7 kpc
estimated from the formula provided by Allen (1973), where the
relation between the extinction of the star light near the galactic
plane and the distance is estimated based on the average properties of
interstellar absorbing clouds and the grains between them.  We
obtained similar distance by comparing directly the $n_H$ value of
X 1908+075 with that of another galactic plane source GRS 1915+105, which is
$3^\circ$ away from X 1908+075 and is believed to have a kinetic distance of 12.5 kpc
and $n_H = 5\times 10^{22}$ cm$^2$ (\markcite{chaty96}Chaty et al. 1996).  At the distance of 7 kpc,
the estimated X-ray luminosity of X 1908+075 is $L_{1.5-100 keV} \sim
5 \times 10^{36}$ erg s$^{-1}$.

Finally, the detected $180$-d period seems to be caused by scattered solar X-ray contamination.  We found from the ASM light curve that three major flares happened at the times when the sun was near the field of view of this source, while flares with smaller amplitudes are evident when the sun is  $180^\circ$ away from the source position. The latter may be caused by the back-scattering of the solar X-rays off  the earth's atmosphere and the instrument.

\section{Conclusions}

A  persistent  hard X-ray emission ($\sim 8$ mCrab at 1.5--12 keV) of X 1908+075 has been detected with  the {\it RXTE}/ASM at the Einstein/IPC position of $ ra=19^h 10^m 46^s$, \ \ $dec=7^\circ 36'
07''$ (J$2000)$ for a 3-year  period. A  $4.4$-d period has been dicovered from the ASM light curves.  The modulation is nearly sinusoidal, persistent and coherent for  3 years. The amplitude of the best-fit sinusoidal wave is $26 \pm 3 \%$ in the $1.5$--$12$ keV band.  Our results support the simple interpretation that X 1908+075 is  the same source as 4U 1909+07,  3A 1907+074, 1H 1907+074, GPS 1908+075, 1E 1908.4+0730, as well as the source detected by HEAO A-3.  We argue that the 4.4-d period is the orbital period of a high mass binary system.  We also suggest that this system consists of a compact object and a supergiant star that produces a strong stellar wind. The X-ray modulation is then caused by changes in the optical
depth along the line of sight to the X-ray source as a function of orbital phase.

\acknowledgments
We are very grateful to the entire $\it{RXTE}$ team at MIT for their support. Support for this work was provided in part by NASA Contract NAS5-30612.


\section{Introduction}
Lattice\cite{wil} formulation is a gauge invariant regulator of non abelian
gauge
theories. Numerical simulations on the lattice produce from first principles
regulated correlators of physical quantities.

Simulations can be used to compute quantities involving low energy modes, which
are out of reach of perturbation theory, such as weak interaction matrix
elements, masses, matrix elements of operators in the light cone expansion. The
typical problems encountered in this ``phenomenological'' use of lattice are the
removal of the cut off (renormalization), and limitations in computer power.

Simulations can also be used to test theoretical ideas and to investigate the
structure of the theory. An example of investigation of ``theoretical'' type is
the study of the mechanism of confinement. The typical difficulty in this
approach is to have good theoretical ideas to
test numerically, and possibly to falsify from the first principles.

Apart from confinement itself there are a few fundamental issues at the
background of our understanding of QCD. Among them
\begin{itemize}
\item[1)]
The $1/N_c\to 0 $ limit. The conjecture\cite{thooft} is that the
basic properties of a gauge theory, e.g. confinement, are already
contained in the limit in which the number of colours $N_c$ goes
large, with $g^2 N_c$ fixed. Corrections ${\cal O}(1/N_c)$ can be
treated as a small perturbation. A consequence of this conjecture
is that also quark loops can be viewed as a small perturbation,
apart from their effect on the scale of the theory. Indeed apart
from the loop with two vertices which is proportional to $g^2
N_f\sim N_f/N_c\sim 1$, and which enters in the $\beta$ function,
loops with more vertices have additional factors $1/N_c$ and are
negligible. According to this conjecture also the mechanism of
confinement and the corresponding order parameter should then be
marginally affected by the presence of quarks.
\item[2)] Understanding the ground state is also important to understand why
perturbation theory works at small distances.
Perturbative quantization describes interaction of quarks and gluons, and the
ground state is the Fock vacuum. Quarks and gluons are not observed in nature,
and the Fock vacuum is certainly not the ground state. This reflects in the
renormalized perturbation expansion as a lack of convergence, even in the sense
of asymptotic expansion\cite{mul}.
\end{itemize}
\section{Confinement in Nature.}
Colour is confined in nature. The expected ratio
of abundance of quarks $n_q$ to abundance
of nucleons $n_p$ is in the standard cosmological model\cite{Okun}
\begin{equation}
\frac{n_q}{n_p} \simeq 10^{-12} \ . \label{eq:e1}\end{equation}
The experimental upper limit is
\begin{equation}
\frac{n_q}{n_p}\leq 10^{-27} \ , \label{eq:e2}\end{equation}
coming from Millikan like experiments on $\sim 1\,g$ of matter.

The estimate (\ref{eq:e1}) is conservative. If we assume no confinement and $T$
is the temperature at which quarks decouple their effective mass $m_q$ is
$\sim T$. The reactions
\[ q + \bar q \to {\rm mesons}\qquad
q + q \to \bar q +\,{\rm baryons}\]
are esothermic: let $\sigma$ the corresponding cross section and
$\sigma_0 \equiv \lim_{v\to 0} \sigma v$. Then quarks will decouple when
\begin{equation}
n_q \sigma_0 = G_N^{1/2} T^2 \ . \label{eq:e3}\end{equation}
Since $n_\gamma \sim T^3$ this implies
\begin{equation}
\frac{n_q}{n_\gamma} \sim \frac{ G_N^{1/2}}{T \sigma_0} =
\frac{10^{-18}}{m_p m_q \sigma_0} \ , \label{eq:e4}\end{equation}
$\sigma_0\sim m_\pi^{-2}$. The ratio (\ref{eq:e1}) corresponds to $m_q = T =
10\,GeV$.

The factor $10^{-15}$ between the observation and the expectation cannot be
explained by fine tuning of a small parameter. Like the experimental limit on
the resistivity of a superconductor, it can only be explained in terms of
symmetry.

A suggestive idea in that direction is that vacuum is a dual superconductor
\cite{Mandelstam}.
The chromoelectric field between a $q$ $\bar q$ pair is constrained by dual
Meissner effect into an Abrikosov flux tube with energy proportional to the
distance.

A relativistic version of the free energy of a superconductor, which is the
analog of effective action in field theory, is
\begin{equation}
G = -\frac{1}{4} F_{\mu\nu} F_{\mu\nu} +
(D_\mu\varphi)^* (D_\mu\varphi) + V(\varphi) \ , \label{eq:e5}\end{equation}
where
\begin{equation}
D_\mu\varphi = (\partial_\mu - i q A_\mu)\varphi\label{eq:e6}\end{equation}
is the covariant derivative and
\begin{equation}
V(\varphi) = \mu^2 \varphi^*\varphi - \frac{\lambda}{2}(\varphi^*\varphi)^2
\label{eq:e7}\end{equation}
is the effective potential.
$\mu^2$ and $\lambda$ are funtions of the temperature, and $\mu^2(T) > 0$ in
the superconducting phase, where the potential has a mexican hat shape.

Putting $\varphi = f e^{i q \theta}$, with $f > 0$ gives
$
D_\mu \varphi = - i q (A_\mu - \partial_\mu\theta) f e^{i q \theta}
$.
Under a gauge transformation $A_\mu \to A_\mu + \partial_\mu \Lambda$
\begin{equation}
\partial_\mu \theta \to \partial_\mu \theta + \partial_\mu\Lambda
\label{eq:e9}\end{equation}
and $\tilde A_\mu = A_\mu - \partial_\mu\theta$ is gauge invariant.
Moreover
\begin{equation}
\tilde F_{\mu\nu} = \partial_\mu\tilde A_\nu - \partial_\nu\tilde A_\mu=
F_{\mu\nu}\label{eq:e10}\end{equation}
and the free energy can be rewritten as
\begin{equation}
G = -\frac{1}{4}\tilde F_{\mu\nu}\tilde F_{\mu\nu} +
q^2 f^2 \tilde A_\mu\tilde A_\mu +
\partial_\mu f \partial_\mu f + \mu^2 f^2 - \frac{\lambda}{2} f^4
\label{eq:e11} \ . \end{equation}
The equations of motion are
\begin{equation}
\partial_\mu \tilde F_{\mu\nu} + q^2 f^2 \tilde A_\nu = 0 \ , \qquad
(\partial^2 + \mu^2) f = 2 \lambda f^3
\label{eq:e12}\end{equation}
A static solution in the gauge $A_0=0$ has $\vec E = 0$, $f = \tilde f \equiv
\sqrt{\mu^2/2\lambda}$ and eq.(\ref{eq:e12}) reads
\begin{eqnarray}
\vec\nabla\wedge\vec H + q^2 {\tilde f}^2 \vec{\tilde A} &=& 0
\label{eq:e13a} \ , \\
\nabla^2\vec H - q^2 {\tilde f}^2 \vec H &=& 0  \ . \label{eq:e13b}
\end{eqnarray}
Eq.~(\ref{eq:e13a}) means that also in the absence of electric field there is a
permanent current, or that $\sigma = \infty$. Eq.~(\ref{eq:e13b}) is Meissner
effect: the field $\vec{\tilde A}$ penetrates by a length $(q \tilde f)^{-1}$.
At large distance from the center of a flux tube $\vec{\tilde A} = 0$ or
\begin{equation}
e \oint \vec{\tilde A} d\vec x = e \Phi(B) = 2\pi n \ ,
\label{eq:e14}\end{equation}
which is the Dirac quantization condition.
Abrikosov flux tubes have monopoles at their ends.

Those phenomena are a consequence of symmetry\cite{wei}:
the order parameter is $f = \langle | \varphi| \rangle$, or the non vanishing
v.e.v. of a charged
operator. For dual superconductivity the signal of the phase should be the
v.e.v. $\langle\mu\rangle$ of an operator carrying magnetic charge.
\section{Phenomenology of confinement on the lattice.}
Lattice produces evidence for confinement. Wilson loops, defined as parallel
transport along square contours in space time, provide the static force between
$q\bar q$ pairs, in the limit of large $T$
\begin{equation}
W(R,T) \mathop{\simeq}_{T \to \infty}
\exp(-\sigma V(R) T)\label{eq:e15}\end{equation}
The area law observed in lattice gauge theory\cite{creu}
\begin{equation}
W(R,T) \simeq \exp(-\sigma R T)\label{eq:e16}\end{equation}
means
\begin{equation}
V(R) = \sigma R \label{eq:e17}\end{equation}
or that an infinite amount of energy is needed to pull the two particles at
infinite distance from each other. $\sigma$ is the string tension, related
to the slope of Regge trajectories.

Also chromoelectric flux tubes between $q\,\bar q$ pairs are observed, with
transverse size $\sim 0.5\,{\rm fm}$\cite{hay,adg1}. The colour orientation of
the chromoelectric field inside them can also be studied\cite{adg1,gree}.

Finally the collective modes of the string formed by the flux tube can be
analysed\cite{Gliozzi}.

All these facts support the picture of confinement as due to dual
superconductivity of vacuum.
A microscopic understanding is however needed. In particular monopoles
which condense in the vacuum have to be identified.

\section{Monopoles}
In QED, which is a $U(1)$ gauge theory, magnetic charges are omitted,
since they are not observed in nature. As a consequence the general solution
of Maxwell's equations can be given in terms of a vector potential $A_\mu$.
The field strength tensor
\begin{equation}
F_{\mu\nu} = \partial_\mu A_\nu - \partial_\nu A_\mu\label{eq:e18}
\end{equation}
obeys the equations
\begin{equation}
\partial_{\mu} F_{\mu \nu} = j_{\nu} \ .
\end{equation}
The absence of magnetic charges indentically follows from eq.(\ref{eq:e18}).
The dual tensor $F^*_{\mu\nu} \equiv \frac{1}{2}
\varepsilon_{\mu\nu\rho\sigma}F^{\rho\sigma}$ is, by virtue of eq.
(\ref{eq:e18}) identically conserved:
\begin{equation}
\partial_\mu F^*_{\mu\nu} = 0 \ . \label{eq:e19}\end{equation}
Eq.~(\ref{eq:e19}) is known as Bianchi identity.

The only way to have a monopole and to preserve Bianchi identity is
\cite{Dirac} to introduce a singularity, and consider the monopole as the
end point of an infinitely thin solenoid (Dirac string), which can be
made invisible if the parallel transport of any charge $q$ around it
is trivial or if
\begin{eqnarray}
\exp \left(iq \oint \vec{A} \cdot \,d \vec{x}\right) = 1 \ .
\label{eq:e20}
\end{eqnarray}
The line integral is intended on a path which encircles the string and is
equal to the magnetic flux, or to $4 \pi$ the magnetic charge of the monopole.
Eq.~(\ref{eq:e20}) implies $4 \pi q M = 2 n \pi$ or $q M = \frac{n}{2}$, which
is the celebrated Dirac quantization condition. As a consequence the theory
becomes compact.

In non abelian gauge theories, in the familiar multipole expansion the
monopole term obeys abelian equation of motion, has a Dirac string and a
number of independent abelian magnetic charges which is $N-1$ for the gauge
group $SU(N)$\cite{Coleman}.
The 't Hooft-Polyakov\cite{'tHooft,Polyakov} monopole
of the $SO(3)$ Georgi-Glashow model obeys this classification.
\section{Monopoles in QCD}
To understand the monopoles in QCD we shall phrase the classification
of ref.~\cite{Coleman} in the language of ref.~\cite{'tHooft2}, or in terms
of ``abelian projections''. We shall refer to $SU(2)$
for simplicity: extension to $SU(N)$ is trivial\cite{adg2}.

Let $\vec{\phi}(x) \cdot \vec{\sigma}$ be any local operator belonging to the
adjoint representation. We define $\hat{\phi}(x) = \frac{\vec{\phi}(x)}
{\left|\vec{\phi}(x)\right|}$. $\hat{\phi}$ is well defined except at zeros
of $\vec{\phi}(x)$. Consider the field strength tensor\cite{'tHooft}
\begin{eqnarray}
F_{\mu \nu}(x) = \hat{\phi} \cdot \vec{G}_{\mu \nu}(x) - \frac{1}{g}
\hat{\phi}(x) \cdot \left( D_{\mu} \hat{\phi}(x)\wedge D_{\nu}
\hat{\phi}(x) \right)
\label{eq:e21}
\end{eqnarray}
where $\vec{G}_{\mu \nu}(x) = \partial _{\mu} \vec{A} _{\nu} -
\partial _{\nu} \vec{A} _{\mu}$ and $D _{\mu} = \partial _{\mu} - g  \vec{A}
\wedge$ is the covariant derivative. The coefficient of the second
term in eq.~(\ref{eq:e21}) is chosen in such a way that the quadratic term
$\vec{A}_{\mu} \wedge \vec{A}_{\nu}$ cancels with the first term.
Both terms are gauge invariant under regular gauge tranformations.

A gauge transformation $U(x)$ which brings $\hat{\phi}$ along the 3 axis
$\hat{\phi} = (0,0,1)$, and diagonalizes $\vec{\phi}(x) \cdot \vec{\sigma}$
is called an abelian projection. After abelian projection
\begin{eqnarray}
F_{\mu \nu} = \partial _{\mu} A^3_{\nu} -  \partial _{\nu} A^3_{\mu}
\end{eqnarray}
is an abelian field. This holds in all points where $U(x)$ is regular.
Defining the dual tensor $F_{\mu \nu}^{\star}$ as
\begin{eqnarray}
F_{\mu \nu}^{\star} = \frac{1}{2} \epsilon _{\mu \nu \rho \sigma}
F_{\rho \sigma} \ ,
\end{eqnarray}
the magnetic current $j_{\mu}^M$ is defined as
\begin{eqnarray}
\partial_{\mu} F_{\mu \nu}^{\star} = j_{\nu}^M
\end{eqnarray}
and is identically conserved. It is zero except at the singular points of
$U(x)$, where monopoles can appear. Thus
\begin{eqnarray}
\partial_{\mu} j^{M} _{\nu} = 0
\end{eqnarray}
defines a magnetic $U(1)$ symmetry of the theory. It is not a subgroup of
the gauge group since both $F_{\mu \nu}$ and $j_{\mu}^M$ are colour singlet.
If the vacuum is not invariant under that $U(1)$ symmetry, monopoles condense
like the Cooper pairs and there is dual superconductivity.

Notice that
\begin{enumerate}
\item there is a magnetic $U(1)$ symmetry for each operator $\vec{\phi}(x)$
in the adjoint representation.
\item Under the abelian projection $U(x)$, due to singularities, the field
strength tensor $\vec{G}_{\mu \nu}$ acquires a singular component
\begin{eqnarray}
\vec{G}_{\mu \nu} \to U \vec{G}_{\mu \nu} U^{-1} +
\vec{G}_{\mu \nu} ^{sing} \ .
\end{eqnarray}
The regular part can have monopole sources. $\vec{G}_{\mu \nu} ^{sing}$
describes Dirac strings starting from the monopoles.
\end{enumerate}

The strategy will then be to detect condensation of different monopole
species by measuring with numerical simulations across the deconfining
transition a ``disorder'' parameter, which detects dual superconductivity.
The usual order parameter is the Polyakov line. Our disorder parameter
will be zero in the deconfined phase, where the order parameter is
non zero, and different from zero in the confined phase, where it is zero.
The concept of disorder parameter is typical of systems which admit a dual
description\cite{Kadanoff}. They usually have extended structures
with non trivial topology (monopoles in QCD), which condense in the disordered
phase. In a dual description these structures are described by local fields,
and the definition of order and disorder is interchanged. Before giving the
results and a few details on how they have been obtained, we shall present
the expectations.
\section{Expectations vs. results}
As we have seen dual superconductivity of the vacuum is not a well defined
concept. There are infinitely many choices for the operator $\vec{\phi}(x)$,
and for each of them ther can or can not be condensation. What is the good
choice, if any?
\begin{enumerate}
\renewcommand{\theenumi}{\Alph{enumi}}
\item There is club of practitioners of the ``maximal abelian projection'',
saying that their choice is better than others. In fact with this choice the
abelian field ``dominates'' the configurations, in particular the part of
it which is produced by monopoles. This can prove convenient to attempt a
construction of effective lagrangeans, but in principle does not preclude
any pattern of symmetry.
\item There is a conjecture that all abelian projections are
equivalent\cite{'tHooft}.
\end{enumerate}
Most probably both attitudes reflect our imperfect knowledge of the symmetry
of the disordered phase.

Discriminating between (A) and (B) is possible on the lattice. The results
that we have obtained by a systematic study of dual superconductivity
show unambiguously
\begin{enumerate}
\item That confinement is a transition from normal to dual superconductor
ground state. This is strong evidence that the mechanism of confinement is
indeed dual superconductivity of the vacuum.
\item A few different abelian projections have been analyzed. All of them
show the same behaviour. The scenario (B) seems to be true.
\end{enumerate}

This is reassuring for the validity of the mechanism itself.
If only one abelian projection would show dual superconductivity
only the particles with non zero charge with respect to that
$U(1)$ could be confined:
there exist states for which that charge is zero, e.g. the gluon which is
parallel to $\hat{\phi}$. Moreover the colour direction of the
electric field in the flux tubes observed in the lattice should also be
parallel to $\hat{\phi}$, being the electric partner of the magnetic $U(1)$
of the monopoles. This has been shown to be not true\cite{adg1,gree}.
Both these facts are naturally explained if the scenario (B) is at work.

We conclude by giving a few details on the technique used to detect dual
superconductivity. The technique has been checked in many well known
systems showing order disorder duality versus traditional
descriptions\cite{tutti}.
\section{The disorder parameter}
An operator $\mu$ is constructed which carries non zero magnetic charge
with respect to the $U(1)$ under study.
Finite temperature is realized
on the lattice by the usual thermodynamical recipe of having euclidean
time running from $0$ to $1/T$ (the temperature), with periodic boundary
conditions for barions, antiperiodic for fermions. On a lattice this is
done by using a size $N_s^3 \times N_t$, with $N_s \gg N_t$ and $a(\beta) N_t
= 1/T$ ($\beta = 2N/g^2$). By renormalization group arguments
$a(\beta) \simeq \frac{1}{\Lambda _L} \exp(- b_0 \beta)$,
or $T = \frac{\Lambda _L}{N_t} \exp( b_0 \beta)$.

The technique used to construct $\mu$ is inspired to ref.~\cite{Kadanoff}
and~\cite{Swieca} and is a complicated version of the simple formula for
translations in elementary quantum mechanics
\begin{eqnarray}
e^{ipa} | q \rangle = | q + a \rangle \ .
\end{eqnarray}
In the Schr\"odinger representation, the field $\vec{A}(x,t)$ plays the r\^ole
of $x$, the conjugate momentum $\vec{\Pi}(\vec{x},t)$ the r\^ole of $p$ and
\begin{eqnarray}
\mu |\vec{A}(\vec{x},t) \rangle \equiv \exp^{\left( i \int \,d^3 x \vec{\Pi}
(\vec{x},t) \cdot \vec{\bar{A}}( \vec{x}- \vec{y}) \right)}
| \vec{A}(\vec{x},t) \rangle = | \vec{A}(\vec{x},t) + \vec{\bar{A}}(\vec{x} -
\vec{y}) \rangle \ .
\label{eq:e31}
\end{eqnarray}
If $\vec{\bar{A}}$ is the field of a monopole at $\vec{y}$, $\mu$ is indeed
the creation operator for a monopole, at site $\vec{y}$ and at time $t$. When
inserted in the Feynman integral the operator is nothing but a linear term
in the conjugated momentum added to the lagrangean and hence a shift of $
\vec{\Pi}$ at
time $t$. Care is needed to adapt the definition (\ref{eq:e31}) to a compact
system, in a form which does not depend on the choice of the gauge for the
classical field $\vec{\bar{A}}$. This can be done, and the result is, as
sketched above\cite{adg2,tutti}
\begin{eqnarray}
\langle \mu \rangle = \frac{Z[S + \Delta S]}{Z[S]} \ ,
\end{eqnarray}
with $\Delta S$ different from zero on a hyperplane at constant $x_0$.

Being the exponential of a sum on $N_s^3$ sites $\langle \mu \rangle$ is
subject to strong fluctuations. Numerically it is better to measure
the quantity
\begin{eqnarray}
\rho = \frac{\mbox{d}}{\mbox{d}\beta} \log \langle \mu \rangle =
\langle S \rangle _S -  \langle S + \Delta S \rangle _{S + \Delta S}
\end{eqnarray}
and to reconstruct $\langle \mu \rangle$ as
\begin{eqnarray}
\langle \mu \rangle = \exp \left[ \int_0^{\beta}
\rho(\beta^{\prime}) \,d \beta^{\prime} \right] \ .
\end{eqnarray}
A typical shape of $\rho$ is shown in fig.~\ref{fig1}, one of $\langle
\mu \rangle$ in fig.~\ref{fig2}. As is well known\cite{Lee} $\langle \mu
\rangle$ as an analytic function of $\beta$, can not vanish indentically in
the deconfined phase if the number of degrees of freedom is finite.
An extrapolation to $N_s \to \infty$ must be done by finite size analysis.
Fig.~\ref{fig3} shows that a few different abelian projections behave
in the same way. Fig.~\ref{fig4} shows that two different monopole
species of $SU(3)$ also behave in the same way.

\noindent
\begin{figure}
\begin{minipage}[htb]{7.8cm}
\centering{\epsfig{figure=fig2.eps, width=7.8cm}}
\caption{$\rho$ vs. $\beta$ for $SU(2)$ gauge theory. Plaquette projection,
lattice $16^3 \times 4$.}
\label{fig1}
\end{minipage}
\hfill
\begin{minipage}[htb]{7.8cm}
\centering{\epsfig{figure=fig1.eps, width=7.8cm}}
\caption{$\langle \mu \rangle$ vs. $\beta$ for $SU(2)$ gauge theory.
Plaquette projection, lattice $16^3 \times 4$.}
\label{fig2}
\end{minipage}
\end{figure}

As for the extrapolation to $N_s \to \infty$ three different ranges
of $\beta$ are considered:
\begin{enumerate}
\item $\beta \gg \beta_c$: there perturbative theory is at work and $\rho
\simeq - c N_s + d$, with $c$ a positive constant. As $N_s \to \infty$
$\rho  \to - \infty$ and $\langle \mu \rangle = 0$.
\item $\beta \ll \beta_c$: there $\rho$ has a finite limit as $N_s \to \infty$,
and $\langle \mu \rangle$ is finite.
\item $\beta \simeq \beta_c$: there we expect $\langle \mu \rangle \simeq
\left( \beta_c - \beta \right) ^{\delta}$. From dimensional analysis
\begin{eqnarray}
\langle \mu \rangle = N_s^{\delta} \Phi \left( \frac{a}{\xi}, \
\frac{\xi}{N_s}, \ \frac{N_t}{N_s} \right) \ .
\end{eqnarray}
Where $\xi \gg a$, and if $N_s \gg N_t$
\begin{eqnarray}
\langle \mu \rangle = N_s^{\delta} \Phi \left( 0, \
\frac{\xi}{N_s}, \ 0 \right) \ ,
\end{eqnarray}
and since $\xi \simeq \left( \beta_c - \beta \right)^{- \nu}$, we get
the scaling law
\begin{eqnarray}
\frac{\rho}{N_s^{1/\nu}} = f \left( N_s^{1/\nu}(\beta_c - \beta)\right) \ .
\end{eqnarray}
The quality of the scaling is shown in fig.~\ref{fig5} for $SU(2)$ and
in fig.~\ref{fig6} for $SU(3)$. It gives a determination of
the critical indices $\nu$, $\delta$, and of $\beta_c$.
For $SU(2)$ we get
\begin{eqnarray}
\nonumber
\hspace{4cm}
\begin{array}{l}
\nu = 0.63(5)\\
\beta_c = 2.30(2) \qquad N_t = 4 \\
\delta = 0.20(8)
\end{array}
\ .
\end{eqnarray}
$\nu$ and $\beta _c$ are in agreement within errors  with independent
determinations, and $\nu$ indicates a second order phase transition.

For $SU(3)$
\begin{eqnarray}
\nonumber
\hspace{4cm}
\begin{array}{l}
\nu = 0.33(2)\\
\beta_c = 5.70(3) \qquad N_t = 4 \\
\delta = 0.54(4)
\end{array}
\ ,
\end{eqnarray}
indicating that the transition is first order.
\end{enumerate}

\noindent
\begin{figure}
\begin{minipage}[htb]{7.6cm}
\centering{\epsfig{figure=fig3.eps, width=7.6cm}}
\caption{$\rho$ vs. $\beta$ for different abelian projections. $SU(2)$ gauge
theory, lattice $12^3 \times 4$.}
\label{fig3}
\end{minipage}
\hspace{0.3cm}
\begin{minipage}[hb]{7.5cm}
~\vspace{1mm}
\centering{\epsfig{figure=fig4.eps, width=7.5cm}}
\caption{$\rho$ vs. $\beta$ for different abelian generators of $SU(3)$.
Polyakov projection, lattice $12^3 \times 4$.}
\label{fig4}
\end{minipage}
\end{figure}
\noindent
\begin{figure}
\begin{minipage}[htb]{7.8cm}
\centering{\epsfig{figure=fig5.eps, width=8cm}}
\caption{Quality of the scaling for $SU(2)$ gauge theory.
Plaquette projection, $N_t = 4$.}
\label{fig5}
\end{minipage}
\hfill
\begin{minipage}[htb]{7.8cm}
\centering{\epsfig{figure=fig6.eps, width=7.8cm}}
\vspace{0.2cm}
\caption{Quality of the scaling for $SU(3)$ gauge theory.
Polyakov projection, $N_t = 4$.}
\label{fig6}
\end{minipage}
\end{figure}

The method used has been tested on a number of known systems and understood
in its details. The result show beyond any reasonnable doubt that dual
superconductivity occurs, in different abelian projections, in connection
with confinement.

\vspace{0.5cm}
The part of the work reported due to our group has been done in collaboration
with L. Del Debbio, G. Paffuti, P. Pieri, B. Lucini, D. Martelli in the last
few years. Their contribution was determinant to the results.


\section{Basic Equations, Symmetries}

\noi
When pair creation can be neglected, a system of two scalar particles can be
 described by a pair of coupled Klein-Gordon equations
$$ \displaystyle
2 H_a   \Psi = m _a ^ 2  \  \Psi  \qquad \quad  a,b = 1,2     $$
referred to as the {\it mass-shell constraints\/}.
 Here $\Psi $ has two arguments $q_1, q_2$ running in spacetime.
We  cover all cases of practical interest assuming that
$ H _a = K_a + V $.

\noi        In the above  formula  $2K_a =  (p_a - e_a  A (a) ) ^2$
 is the squared-mass operator for particle $a$ alone in the magnetic field,
 and $V$ is a suitable modification of the term  $ V \zer $
 which would describe the mutual interaction in the absence of external field;
 this modification is necessary in order to
keep the mass-shell constraints mutually compatible when the  field
$F_{\alp \beta} $  is applied.

\noi         For all vectors
        $\xi , \eta $ we write
   $\xi \cdot    F  \cdot \eta $ for
   $\xi ^\alpha     F _{\alpha \beta}   \eta ^\beta  $. With a similar
convention
 $  \disp   A(a) = \half  \     q_a  \cdot  F $
in a  Lorentz-covariant gauge.

\noi   Notice that $\xi \cdot \eta \longi =  \xi \longi \cdot \eta \longi $
and    $\xi \cdot \eta \st  =  \xi \st  \cdot \eta \st  $.

    \noi
An important technical point is that  applying a constant magnetic field
 provides a unique and invariant decomposition of any four-vector $\xi$ into
 longitudinal and transverse parts, say
  $ \xi = \xi \longi  + \xi \st $.   The orthocomplement of $(E\longi)$ in the
space of four-vectors  is a  two-dimensional  plane $(E\st)$ endowed with
elliptic metric. In any adapted frame, $\xi \longi$ (resp. $\xi \st$)
has nonvanishing coordinates $\xi ^0 , \xi ^3 $ (resp. $ \xi^1 , \xi^2$).

\medskip

\noi  The theory of relativistic two-body systems,  formulated  many  years
 ago along the lines of "predictive mechnics" and "constraints theory"
\cite{rep} \cite{tod} \cite{barcelo} \cite{coupl}
has been more recently extended to  cases where some  external field is present
\cite {bij89} \cite{DVNCim}. Here we
assume that a constant homogeneous magnetic field
 is applied to a pair of opposite charges, say       $e_1 =-e_2 = e$.

\medskip
\noi        The constraint approach employed here has over the
 Bethe-Salpeter  equation
several advantages; for example in the particular case of isolated systems
 (no field applied)   the dependence on relative
time  gets automatically factorized out  \cite{sazdj}.

\noi It is convenient to re-arrange the canonical variables as follows
$$ z= q_1 -q_2    \qquad \qquad               Q= \half (q_1 + q_2)    $$
$$ y = \half  \   (p_1 -p_2) \qquad  \qquad     P = p_1 + p_2         $$
so $ [z,y]= [Q,P]=-i \delta $, etc.

\noi  The Lie algebra of the Lorentz group is generated by the tensor
$$      M   =    M_1  + M_2     =     Q\wedge P + z\wedge y      $$
with  $M_1    = q_1\wedge p_1 , \quad   M_2 =   q_2 \wedge p_2$.
In any adapted frame, rotations in $(E \st)$ are   generated by $M_{12}$
and boosts   in $(E \longi )$ by $M_{03}$.

\noi
An essential ingredient of mutual interactions  \cite{rep} is the quantity
$\disp   {\widetilde z}^2 = z^2 - (z \cdot P )^2 / P^2 $.
But in order to avoid denominators in calculations,
it is  convenient to  employ
  \begin{equation}  Z = z^ 2  P ^ 2  -  ( z \cdot  P )^ 2    \label{2.1}  \end{equation}
We shall assume   that
      \begin{equation}  V\zer = f(Z, P^ 2, y \cdot  P)   \label{2.2}   \end{equation}
This form is general enough to accomodate a large class of interactions.

\medskip
\noi  {\sl Definition}

\noi
When speaking of {\em energy-dependent interactions}, we refer to the {\em
total} energy of isolated systems, namely $\sqrt {P^2}$.

\noi               Although $Z$ is more
practical for calculations, it would be more  natural to take
 ${\widetilde z}^2$ and $P^2$ as  independent dynamical variables, defining
$\disp   g ( Z/P^2, P^2 , y \cdot P)   =   f(Z, P^ 2, y \cdot  P) $.
Therefore we  say  that $V\zer$ does{\em not} depend on (total) energy when
the function $f$ takes on the form
$f = g(Z/P^2, y \cdot P) $.

\medskip   \noi
Although  $f$ in (\ref{2.2}) is supposed to be known, it would be a problem
to determine $V$ in closed form.
In the  external-field representation, which involves  a new wave function
 $\Psi '$ and new operators $H_a ',  K_b ', V' $, this problem is solved by
making the ansatz
\begin{equation}    V^\prime  =  f( \widehat Z, P^2, y\longi .P \longi )
     \label{2.3}  \end{equation}
where
$ \widehat Z= Z' \zero=    (Z')_   {F=0}     $ (it turns out that
$     \widehat Z $ commutes with $  y\longi \cdot P \longi $).
The explicit form of $ \widehat Z $ was calculated in ref. \cite{DVNCim}.
      \begin{equation} \widehat Z =
Z  +  2(P\longi ^2 \   z \cdot P - P^2  \   z\longi \cdot P\longi  ) L
  +  P\st^2 P\longi ^2  L^2                    \label{2.4}  \end{equation}
   where the scalar  $L$ is defined as
 \begin{equation}  L = {P \longi \cdot z   \over  (P \longi )^ 2}   \label{2.5} \end{equation}
The equations of motion are compatible provided  \cite{DVNCim} that
 $ \widehat Z $ commutes with $y\longi \cdot P \longi $.
 Let us transform (\ref{2.4}) in order to render this  commutation
property manifest.
First we split  $z$ as the sum of $z\longi $ and $z \st $ in   $Z$, hence
\begin{equation} Z = (z\st ^ 2 + z \longi ^ 2 ) P ^2
- (z \st \cdot P )^ 2 - ( z\longi \cdot P )^ 2
-2 (z \st \cdot P )(z \longi \cdot P )            \label{2.6}      \end{equation}
Develop (\ref{2.4}) and perform  elementary manipulations
 using (\ref{2.6}). We get
$$ \widehat Z = Z + 2 (z \st \cdot P )(z \longi \cdot P)
-(z \longi \cdot P)^ 2   {P \st ^ 2 \over P \longi ^ 2}            $$
Using (\ref{2.6}) we notice cancellation of the terms proportional to
 $     (z \st \cdot P )  (z \longi \cdot P)$ and we  can write
  \begin{equation} \widehat Z =
 z \st ^ 2 P ^ 2 - (z \st \cdot P)^ 2  +
 P^ 2 (  z \longi ^ 2   -   {  (z \longi \cdot P \longi ) ^ 2 \over P
 \longi ^ 2}  )
     \label{2.7} \end{equation}
It is convenient to define  the projector  "orthogonal" to $P \longi$, say
\begin{equation} \Omega ^ \alpha _ \beta =   \delta ^ \alpha _\beta  - { P \longi ^ \alpha
 { P \longi }_\beta   \over P\longi ^ 2 }     \label{2.8}  \end{equation}
   because we can write
\begin{equation}  z \longi ^ 2    -  { (z \longi \cdot P \longi ) ^ 2
 \over P \longi ^ 2}   =  (\Omega z \longi )^ 2              \label{2.9} \end{equation}
and we easily check that  $ (\Omega  z) ^\alpha $  commutes with
    $ (y \longi \cdot P\longi  )$.
So we  have
\begin{equation} \widehat Z = z \st ^ 2  P^ 2 - (z\st \cdot P)^ 2
   +  (\Omega z \longi )^ 2   P^ 2                     \label{2.10}    \end{equation}
which justifies the claim that   $ \widehat Z$
commutes with   $y\longi \cdot P\longi$
Here we notice that
     $ \Omega z \st  = z \st $ and
\begin{equation}   (\Omega z)^2  = z \st ^2  + (\Omega z \longi)^2    \label {2.10bis}
                                                                         \end{equation}
thus  we  finally obtain
\begin{equation} \widehat Z   = (\Omega z)^ 2 P^ 2 - (z \st \cdot P )^ 2
  \label{2.11} \end{equation}
which is much more tractable than formula (\ref {2.4}).

\medskip
\noi
Mass-shell constraints can be replaced by their sum and difference, so we set
$$   \mu = \half (m_1 ^2 + m_2 ^2) ,  \qquad \quad
 \nu = \half (m_1 ^2 - m_2 ^2)                      $$
  The  explicit form of   $K'_1$ and  $  K'_2$
       was given in Ref. \cite{DVNCim}.
  Equations  (3.36) of Ref. \cite{DVNCim}  yield in the present notation
\cite{novecL}
  \begin{equation} K'_1 + K'_2 =
K_1 + K_2  -2T  \   {y \longi \cdot P \longi  \over  P \longi  ^2 }
  + {T^2  \over  ( P \longi ) ^ 2}     \label{2.12}  \end{equation}
where
 \begin{equation} T = K_1 -K_2 -  y\longi \cdot P\longi              \label{2.13} \end{equation}
and the  difference is
\begin{equation}  K'_1 - K'_2 = y\longi \cdot P \longi               \label{difK} \end{equation}
It is noteworthy that  $M_{03} $ and $M_{12}$ are not affected by going over to
the external-field representation. In other words we can write
\begin{equation}     M'_{12} =  M_{12} ,  \qquad    M' _{03} =   M _ {03}
                                                \label {rotboost}     \end{equation}
Indeed the transformation from $\Psi$ to $\Psi'$ is formally generated by
$B= LT$ where $L$ and $T$ are given by  (\ref{2.5})     and     (\ref{2.13})
respectively \cite{rigo}.
  Commutation of $L$ with $ M_{12}$ and   $  M _ {03}$ is
obvious.  For commutation of $T$, the only point to be checked is that
 $K_1 - K_2$    actually    commutes with $M_{12}$. But  $K_a = K (a) $
where  $K(a)$ is the (half-squared) squared-mass operator  for    particle $a$
alone in the field. We know the constants of the motion in the one-body
sector \cite {bacry}.   In particular we know that $K_a$ commutes with
both  $ (M_a)_{03}$    and      $ (M_a)_{12} $ .   Thus $T$
commutes with     $  M_ {03}$    and      $ M_{12}$.
Finally $B$ shares the same property, which formally proves (\ref{rotboost}).
Let us prove the following

\noi {\sl Proposition

\noi   Angular momentum in    $(E\st)$   and   boost     in $(E \longi)$
 are constants of the motion\/}.

\noi   In other words we claim that our  squared-mass operators  both
commute with  the transverse and longitudinal components of the total
angular momentum.

\noi   Working in the external-field representation, all we need is
to prove  that      $  M_ {03}$    and      $ M_{12}$
commute with    both $ K'_a + V'$, or equivalently with
$K'_1 + K' _2 + 2 V'$ and with   $y \longi \cdot P \longi $.
Commutation with  $K'_1$ and $K'_2$ separately is ensured from
  the properties of single-particle motion in the field.
Moreover     $y \longi \cdot P \longi $ is invariant under any spacetime
rotation. The last point to check is whether    $  M_ {03}$    and    $ M_{12}$
actually commute with $V'$. It is sufficient that they commute with all
arguments of $f$ in formula  (\ref{2.3}), which is true because  these three
arguments  are manifestly Lorentz invariant.

\noi For completeness, it is in order to remind here that pseudomomentum,
originally represented by
$$   C= P +  {e\over 2} \  z \cdot F $$
keeps the same expression in the external-field representation ($C' = C$),
and  is also conserved \cite{DVNCim}.

\subsection {Ultimate Representation}
For neutral systems, a further transformation  inspired by the work of Grotch
and Hegstrom \cite{gro}, and similar to a gauge transformation,
                               permits to get rid of the $Q$ variables.
Transforming the  wave function  yields
$ \Psi '' = ( \exp i\Gamma) \Psi ' $
with the help of  the unitary transformation generated by
  \begin{equation}   \Gamma  =  {e \over 2} (z.F.Q)         \label{2.14}        \end{equation}
We  set
  \begin{equation} {\cal O} ^\sharp  = \exp (i \Gamma) \ {\cal O}  \   \exp (-i \Gamma)
 \    \qquad      {\cal O} '' =
 ({\cal O}' ) ^ \sharp    \qquad  \forall  {\cal O}         \label{2.15}   \end{equation}

\noi       The new equations of motion
\begin{equation} (H''_1 + H''_2) \   \Psi ''  =  \mu \   \Psi ''    \label{2.16} \end{equation}
\begin{equation}   y \longi \cdot  P \longi    \Psi '' =
  \nu \   \Psi  ''   \label{2.17}   \end{equation}
may    "look like"  translation invariant, although  {\it they are not}.
The reason is that  pseudomomentum is transformed to $P ^ \alpha $
by (\ref{2.15}),
that is  $C'' =  P$. Of course $P$ is not any more the generator of spacetime
translations.  These transformations now have a generator $P''$ which differs
from $P$ because $\Gamma$ in  (\ref{2.14}) is not translation invariant.
In the ultimate representation considered  here $C''$  generates the
relativistic analog of the so-called "twisted translations" invoked in
\cite{avr}.

\bigskip
\noi   From now on
 {\it we demand that  pseudomomentum be diagonal with a
timelike four-vector   $k^ \alpha $ as eigenvalue}. Instead of
 $C^\alp  \  \Psi =  k^\alp \  \Psi$  we are using our ultimate
representation   and    write
$  P^\alp   \Psi ''= k ^\alp  \Psi ''$.
Combining this requirement with (\ref{2.17})   we obtain
\begin{equation}  \Psi'' = \exp (i k \cdot Q)
\ \exp (i \nu { z\longi \cdot k \longi \over  | k \longi | } )\
 \phi             \label{2.18}   \end{equation}
where  $ \phi  $ depends   on $z$,  but only through its
projection orthogonal to $k \longi $, and additionally depends on $k$ and on
$\nu$ as parameters.
In other words $ \phi = \phi ( \nu ,k,  \varpi  z) $
with     the following notation.

\noi  {\sl Notation}

For all four-vector $\xi$, we define
$\varpi \xi  $ as  the  projection  of  $\xi$ onto  the 3-plane orthogonal to
$k \longi $,  say
$ \disp (\varpi \xi)  ^\alp =
\xi ^\alp -  ( \xi \cdot k \longi) \     k \longi ^\alp   /  k \longi ^2     $.

\noi
Similarly   $\xi  \per$ denotes the  projection of  $\xi$ onto the 3-plane
 orthogonal to  $k$. In general  $\varpi z    \not= z \per$, but they
coincide when  $k \st $ vanishes.

\medskip
\noi It is convenient  to  introduce here the {\it motional parameter\/}
 $ \disp   \epsilon =  {|k\st | \over |k \longi| }$.  When $\epsilon$ doesnot
vanish, a number of terms involving the contraction  $k \cdot F$ arise. In
fact    $ (k \cdot F)^\alp =  |k|  \  E ^\alp $ where   $  E ^\alp $  is the
electric field   "seen" by an inertial observer moving with constant momentum
  $k ^\alp$  ({\it   motional electric field}).
We have  the identity
  \begin{equation}  \disp   {1 \over k \longi ^2} = {1 \over k^2} (1 - \epsilon ^2)
                                                       \label{2.19}     \end{equation}
Notice that  $k \st$ is linear in $\epsilon$ because we can write
 $k\st =  \epsilon \Lambda  k\longi $ where the second rank
tensor $\Lambda$ represent the boost from the direction of $k\longi$ to the
direction of $k\st$ (thus $\Lambda \cdot \Lambda = \delta$).

\subsection{Explicit Formulas}
  The reduced (or internal) wave function $\phi$ must  be determined
 through the "sum  equation"
 $ (H''_1 + H''_2 ) \   \Psi '' =   \mu   \Psi ''  $,
 simplified with help of   (\ref{2.18}).

\noi            Given the function $f$ involved in (\ref{2.2}),
 let us display $H''_1 + H ''_2 $ in detail.
It is clear that
\begin{equation}   H''_1 + H''_2 =
K''_1 + K'' _2 + 2 V''                               \label{3.1}            \end{equation}
so we have to transform    $  (K'_1 + K'_2)$  and  $V'$
according to (\ref{2.15}). We find   that $Q$ and $z$ are unchanged whereas
\begin{equation} P ^ \sharp = P  + {e \over 2}  F \cdot z
\qquad    \quad          P ^\sharp  \longi   =  P \longi
    \label{3.2}  \end{equation}
  \begin{equation} y^ \sharp = y  - {e \over 2}  F  \cdot  Q     \label{3.3} \end{equation}
\begin{equation} P^ {\sharp 2} = P^2 + e P \cdot F \cdot z
+ {e^2 \over 4} (F \cdot z)^2                      \label{3.4}  \end{equation}
\begin{equation} (K_1 + K_2) ^\sharp =  { P^2 \over 4} + y^2
- {e \over 2} z \cdot   F  \cdot P  + {e^2 \over 4}   (z \cdot F) ^2
                                            \label{3.5}  \end{equation}
\begin{equation} T ^ \sharp =
y \st \cdot P \st - 2e z \cdot F \cdot y    \label{3.6}     \end{equation}
               \noi
Now we apply  transformation (\ref{2.15}) to (\ref{2.12}),
taking (\ref{3.5})(\ref{3.6}) into account.  It gives
 $$  K''_1 +  K'' _2 =      $$
\begin{equation}          {P^ 2 \over 4} + y^ 2   -  {e\over 2}  z \cdot F  \cdot P
 + {e^ 2 \over 4}  (z \cdot F )^ 2
+{T^\sharp \over P\longi ^2} (T^\sharp - 2 y\longi  \cdot P \longi)
                                                          \label{3.7}      \end{equation}
with $T^\sharp $ given by (\ref{3.6}).
We know that $2 V''$ must be added to this expression in order to obtain
 $H''_1 + H''_2  $.
  But in (\ref{2.14})  $F ^ {\mu \nu} $ is purely transverse, therefore
   $ (y \longi \cdot P \longi)^ \sharp  = y \cdot P\longi    $.
We have by  (\ref{2.3})
\begin{equation} V'' = f( {\widehat Z} ^ \sharp,
 {P^\sharp }^2,    y \longi \cdot P\longi )   \label{3.8}  \end{equation}
where $ {P^\sharp }^2 $ is as in (\ref{3.4}) and we must compute
 $\widehat Z ^ \sharp $ from (\ref{2.10}) with   help of (\ref{3.2}).
(We make the convention that
  $\widehat Z ^ \sharp   = ( \widehat Z) ^ \sharp $ and not the reverse).

\medskip
\noi  To this end we apply the transformation (\ref{2.15}) to eq. (\ref{2.10}).
 A glance at  (\ref{2.9})
shows that  $ (\Omega z\longi )^ 2 $ is not affected by the transformation.
Remind that $z$ is unchanged;
we notice that  $ z \st \cdot P ^ \sharp =  z \st \cdot P $ because,
 $F$ being purely  transverse,
$ z \st \cdot  F  \cdot z $ identically vanishes.
Thus, using  (\ref{2.10bis}) we obtain
\begin{equation} \widehat Z ^ \sharp =
{P^ \sharp}^ 2  (\Omega z )^ 2- (z \st  \cdot  P )^ 2     \label{3.9}  \end{equation}
Now, eqs (\ref{3.1})(\ref{3.7})(\ref{3.8}) supplemented with
 (\ref{3.4})and (\ref{3.9}) furnish the complete expression of
  $ H''_1 + H''_2 $,  to be inserted into (\ref{2.16}).
At this stage we are in a position to carry out  the reduction.

\bigskip

\noi   \section{Three-Dimensional Reduction}

\subsection{Calculations}
After transformation to the ultimate  representation we have obtained
  $C'' = P$.

\noi Calculations can be organized as follows:
  Whereas (\ref{2.17})
fixes the dependence in the relative time, eq.(\ref{2.18}) allows us
to factorize out the "center-of-mass motion", and we are left with the reduced
wave function    $ \phi$     which arises in eq. (\ref{2.18}).
Obviously  (\ref{2.17}) implies that
\begin{equation} y \longi \cdot k \longi  \  \phi  = \nu  \      \phi
     \label{4.1}  \end{equation}
thus  $\phi$  depends on $z $ only through its projection $ \varpi z$.
 It is clear that $\phi$ generally depends  on $ \nu $ and $ k$ as parmeters.

\noi
In search for a   reduced wave equation, we  replace $P^ \alpha$
 and  $  y \longi  \cdot P \longi $     respectively
 by  their  eigenvalues  $k^ \alpha $     and  $ \nu $
in $H''_1 + H''_2 $,
 and we divide  by exponential factors.

\noi
For any operator ${\cal O} $ it is convenient to use the following convention
\begin{equation}    ({ \cal O})_{\nu , k } =
{\cal O}   |   _{ y\longi \cdot P \longi = \nu , \    P=k }
\label{4.2}  \end{equation}
The subscript $k$ refers to the vector $k$, which finally contributes by its
longitudinal piece only.
 In this procedure, a term like $y^2$  must be  written  as
    $ \displaystyle  y^2 \equiv  (\Omega y  )^2 +
   { ( y \longi  \cdot  P \longi )^2  \over P \longi ^2 }  $.
If we now introduce the projector $\varpi $ orthogonal to $k \longi $
and  use   identity  (\ref{2.19}) we obtain  for instance,
with help of (\ref{4.1})
\begin{equation}   ({P^2 \over 4 } + y ^2 ) _{\nu, k} =
  { k  ^2 \over 4} + (\varpi y)^2  + {\nu ^2  \over k \longi ^2 }
=    { k  ^2 \over 4} + (\varpi y)^2  + {\nu ^2  \over k^2}
 - \epsilon ^2 {\nu ^2 \over k ^2 }                \label{4.3}   \end{equation}
 which is to be taken into account when computing
 $ (K''_1 + K''_2) _ {\nu , k } $ from (\ref{3.7}).

\noi                 According to  (\ref{3.1})   we have
$ \displaystyle    ( H''_1 + H ''_2)  _{ \nu , k  } =
(K''_1 + K'' _2 )_{\nu , k }  +
     2  (V'') _{\nu , k   } $.
Defining
 \begin{equation}  R (\nu, k \longi, k \st) =   (K''_1 + K''_2)_{\nu, k  }
 \label{4.4}    \end{equation}
\begin{equation} W (\nu, k \longi, k \st) =     (V'')_{\nu,  k }      \label{4.5}   \end{equation}
Recalling  (\ref{3.1}),   equation  (\ref{2.16}) gets  reduced to
\begin{equation}     R    \    \phi  +    2W   \  \phi  =  \mu \phi    \label{4.6}  \end{equation}
\noi Let us stress  that $\mu $ is just a parameter  fixed from the outset.
As other parameters arise in (\ref{4.6}), namely $k$ and $\epsilon$,
the question wether (\ref{4.6}) can be considered as a spectral problem,
 {\it and for which eigenvalue\/}, is not yet settled and will be considered
later on, with help of equations (\ref{4.10})(\ref{4.15}).
See eq. (\ref{4.14}) below.

\medskip
\noi  Since $\phi$ depends on  $z$ only through $\varpi z$, it is important to
realize that     neither   $R$ nor  $W$  involve  the operator
 $z\longi  \cdot k \longi $. This will be  checked below and will permit us to
consider equation (\ref{4.6}) as a three-dimensional problem involving
operators $R$ and $W$ acting on functions of $\varpi z$.

\medskip
 \noi    The explicit expression of $R$ comes from (\ref{3.7}), with help of
(\ref{4.4}).  Since $K''_1$ and $K''_2$ are no more than quadratic in the field
strenght,  let us make the convention that the superscripts $ (1), (2)$
respectively refer to  the (homogeneous) linear and quadratic terms in the
field.  We  start from (\ref{3.7}), compute $K''_1 + K''_2$ to be
inserted into (\ref{3.1}) and further simplify  with help of
convention (\ref{4.2}).
The zeroth order contribution to  $R$ is
  \begin{equation} R\zer = {k^2 \over 4} + {\nu ^2 \over k\longi ^2}   + (\varpi y)^2   +
y\st \cdot k\st \   {y\st \cdot k\st - 2 \nu \over  k\longi ^2}
  \label{4.7} \end{equation}
Applying again identity (\ref{2.19})
and  setting
 \begin{equation} (S)_{\nu , k}  = (\varpi y )^2  +
(y \st \cdot k \st) \
  { y \st \cdot k \st - 2  \nu    \over  k \longi ^2 }
- \epsilon ^2  {\nu ^2  \over k^2 }   \label{4.8} \end{equation}
we can write
 \begin{equation} R   \zer  =
  {k^2 \over 4 } +  {  \nu ^2 \over k ^2} + (S)_ {\nu, k}       \label{4.9} \end{equation}
It is convenient to define
 \begin{equation}  \lam = {k^2 \over 4}  + {\nu ^2 \over k ^2} - \mu      \label{4.10}  \end{equation}
so we can write
\begin{equation} R \zer =     \lam + \mu   + (S)_{\nu , k }
          \label{4.11} \end{equation}
The field-depending terms in (\ref{3.7}) provide
   \begin{equation} R^ {(1)}=  4 e (z\cdot F \cdot y)    {\nu  \over  k \longi ^ 2 }
         -  {e \over 2} z \cdot F \cdot k      \label{4.12}     \end{equation}
\begin{equation} R  ^ {(2)} =
 {e^ 2 \over 4} (z \cdot F)^ 2
+  4 e ^ 2   { (z\cdot  F  \cdot  y)^ 2  \over k \longi ^ 2 }
\label{4.13}       \end{equation}
We remember that $F$ is purely transverse. Contractions involving $F$
only depend on the transverse components;  for instance  $F \cdot k $  is
just a combination of the quantities  $k\st ^\alpha$.
It is noteworthy that only the transverse components of $z, y $ arise in
$ \   R  ^{(1)}, \   R ^{(2)}  \   $, whereas $ (S)_ {\nu, k } $  depends on
  $ \varpi y $ and   $ y \st $.
As a whole,  $R$ depend only on  $ \varpi z $ and  $ \varpi y  $ (recall $y
\st , z \st $ are pieces of $\varpi y , \varpi z $ respectively).

\noi
 In view of (\ref{4.11})(\ref{4.12})(\ref{4.13}), equation (\ref{4.6}) may
be finally written
\begin{equation}  \lambda \phi + [ (S)_{\nu, k} + R ^{(1)} + R ^{(2)} + 2 W  ] \phi  = 0
 \label{4.14} \end{equation}
The square bracket in (\ref{4.14}) is nothing but $ (-N'')_{\nu k}$ provided,
in the original representation,  we introduce the conserved quantity
   \begin{equation}    -N   =
{1 \over 4 } C ^2  +    (C^2)^{-1}   (H_1 - H_2)^2
                     -  (H_1 + H_2)       \label{4.15}              \end{equation}
now represented by   the operator
  $$-N'' =
{P ^2  \over 4} +  {(H''_1 - H''_2) ^2   \over P ^2}
  -  (H''_1 + H''_2)                                           $$
and intimately related with the energy of relative motion.

\medskip
\noi
The last term to be evaluated in (\ref{4.14}) is  $W$.
 In view of (\ref{4.5}) we have first to write down the
expression for $V''$, say (\ref{3.8}).
It follows  that
\begin{equation} W= f(   (\widehat Z ^ \sharp )_{\nu, k} , \
     (P^ \sharp )^ 2 _ {\nu , k} ,  \   \nu )
         \label{4.16}    \end{equation}
 In this  formula $ (P^ \sharp )^ 2 $ is given
by (\ref{3.4}) and  $\widehat Z ^ \sharp $ by (\ref{3.9}).
 Making the substitutions
$P \rightarrow  k$ and  $ y \longi \cdot P \longi  \rightarrow  \nu $,
 hence   $\Omega  \rightarrow  \varpi     $,  we obtain
\begin{equation}  (\widehat Z ^ \sharp )_{\nu, k} =
  (P^ \sharp )^ 2 _ {\nu , k} \   ( \varpi z)^ 2 - (z \st \cdot k)^ 2
      \label{4.17}    \end{equation}
\begin{equation}  ( {P^ \sharp}  ^ 2) _ {\nu, k} =   k^ 2 +
e  \    k \cdot F \cdot z    +{e^ 2 \over 4}(F \cdot z)^ 2
  \label{4.18}  \end{equation}
It is clear that $W$ does not involve the operator $z \cdot k \longi$.
Formulas   (\ref{4.17}) (\ref{4.18}) are to be  inserted into  (\ref{4.16}),
then the explicit form of $W$ will come out.

\medskip
It is natural to consider (\ref{4.14}) as an equation for the eigenvalue
 $\lam$.
But we meet a complication because  $\lam$ is not independent from $k^2$.
In fact we can solve (\ref{4.10}) for $k^2$ and insert the result
      \cite{oneroot}                   into $(N'')_{\nu , k}$.
As a result (\ref{4.14}) bears a  nonlinear dependence on $\lam$. A similar
situation was  pointed out by Rizov, Sazdjian and Todorov \cite{rst} in the
 case of isolated systems undergoing energy-dependent interactions.
 In the presence of  magnetic field
however, the reduced wave equation is nonlinear in $\lam$, even in the simple
case where the mutual interaction term $V\zer$ does not depend on $P^2$.
This can be seen as follows:
first we notice that the occurence of  $ (\widehat Z ^ \sharp )_{\nu, k} $ in
$W$ brings out a dependence on $k^2, k\longi ^2$. Second we observe an
unescapable dependence on    $k^2, k\longi ^2$  in formulas
(\ref{4.8})(\ref{4.12})(\ref{4.13}).
We end  up  with a  nonconventional spectral problem which  requires a
special  treatment, reserved for a future work.

  \medskip
  \subsection{Discussion}

\noi          Finally the mass-shell constraints have  been reduced to the
three-dimensional problem of solving (\ref{4.14}). This formula is nonlinear
in the field strenght and might  be applied to strong fields \cite  {pilk}.
 Let us review  the various contributions it contains.
We distinguish motional terms, depending on $\epsilon$ or depending on
$k\st$, where we know that $k\st$ is linear in $\epsilon$.

\noi  Loosely speaking we could say that,
in as much as the shape of  $W$ departs from the original form assumed
by  $V  \zer$,  every thing goes as if the mutual interaction were somehow
 "modified by the presence of  magnetic field".

\medskip
\noi
{\bf  a) system at rest}

\noi     The  particular case where  pseudomomentum is purely longitudinal
(say $k\st =0$) enjoys a particular simplicity.
If we assume for a moment that $k$  coincides with $k  \longi$,
it is possible to find a frame where ${\vec k}$ vanishes whereas the
electromagnetic field is purely magnetic.
   We refer to this situation as the case  {\it at rest}.

\noi
In this case, $\varpi z = z \per , \  \varpi y = y \per$
and $(S)_{\nu , k} $  simply reduces to $y\per ^2$,
since  $ k\longi $   coincides with $k$.

\noi As $z \st \cdot k$ in  (\ref{4.17})  vanishes, we notice that
                $       (\widehat Z ^ \sharp)   /
  {P^ \sharp} ^ 2 ) _{\nu, k}  $
reduces to  $z \per ^2 $.     According to (\ref{4.16})
 and to a notation defined in Section 2, we can write
$$ W =  g (z \per ^2,  ({P^ \sharp} ^ 2 ) _{\nu, k} ,  \nu)  $$
where $k \cdot F  \cdot z $ vanishes in (\ref{4.18}),  so
$  ({P^ \sharp} ^ 2 ) _{\nu, k}  $ reduces to
$\disp   k ^2  +  {e^2 \over 4} (F \cdot z)^2 $.

\noi  If the mutual interaction doesnot
depend on the energy, we end up with
$ W =  g ( z \per ^2, \nu) $.
Thus, for energy-independent interactions, namely
 $V\zer = g(Z/P^2, y \cdot P)$,
$W$ assumes the form $g (z \per ^2 , \nu) $.
 In other words:

\noi
 {\em   At rest, the magnetic field doesnot modify the mutual interaction,
provided this interaction is  not energy-dependent\/}.

\noi
In contrast, if $\   \disp {\partial V \zer  \over \partial P^2} \   $ doesnot
vanish, the shape of $W$ may  substantially depart from that of  $V \zer$
 in strong fields,  owing to the contribution of
$(F \cdot z )^2    $ in    $  ({P^ \sharp} ^ 2 ) _{\nu, k}$.
This correction to $V \zer$ is a genuine "three-body" term in this sense that
 it vanishes if either the mutual interaction or the magnetic field is turned
off (pretending  that the field is generated by a  ficticious
"third body" located at infinity).

\noi   Looking again at  equation    (\ref{4.14}), we see that,  at rest,
 all surviving terms not included in $W$
can easily be identified as covariant generalizations of the usual terms
present in the non-relativistic theory  \cite{avr} \cite{herold},
 except for a piece of $R^{(2)}$
which depends on the relative angular momentum,  see contribution of
  $ \disp  { z \cdot F \cdot y  \over  k \longi ^2 }  $
 in formula (\ref{4.13}).
 This  contribution  remains
small for heavy systems  ($k ^2 >>  F$)  but  might be significant
for light  systems  ($k^2 <<  F$) in a strong magnetic field.

\medskip
\noi   At {\it first order in the field strenght\/} however,
 the  relative motion admits no
correction other than a term proportional to $\nu$ (indeed $F \cdot k $
vanishes). For equal masses this term is zero and  there is no departure
 from the motion of an isolated system.

\medskip
\noi    {\bf   b)  motional case}

\noi
When $k\st$ is nonzero, we reckognize the  motional electric field contained
in $z   \cdot F \cdot k$.
 For energy-dependent potentials, and even in a  weak field, this term
contributes to  $W$  through  (\ref{4.16}).
  But of course, it may be  neglected in case of slow motion in a weak field,
 where both $\epsilon$ and $F$ are considered as first order quantities,
which  entails that  $F \cdot k$  is  a second order quantity.
On the one hand, this can be seen as a stability property of the
neutral two-body system, under application of a constant field.
 But on the other hand, it forces one to go beyond the
weak-field-slow-motion approximation if one wishes to compute significant
corrections to the energy associated with  relative  motion.

\bigskip
            \section{Conclusion}

\noi      The coupled Klein-Gordon equations describing a globally neutral
system  have been  reduced to a three-dimensional  equation
involving truly motional terms and recoil effects in a   covariant fashion.
In this formulation  the particular symmetry associated with a constant
 magnetic field {\it in space-time\/}  is manifestly respected.
After separation of the internal motion, and after factorizing the
dependence on relative time, the surviving number of degrees of freedom is
finally the same as in the nonrelativistic theory.

\medskip
\noi
 We  now have  a clean theoretical basis for the   study  of relativistic
 bound states in a constant magnetic field, the  simplest of all the  cases
 where  an external field is present.

\noi
In the reduction procedure it was essential to consider eigenstates of
pseudomomentum.
The square of this vector plays the role of an effective squared mass
which can be, in principle, evaluated by solving the reduced wave equation.
But the eigenvalue problem involved in this equation  is crucially
  non-conventional, for the eigenvalue arises in a nonlinear way,
 even if  mutual interaction doesnot depend on the total energy.
This situation  requires a refinement of conventional methods; the method
devised  in Ref. \cite{rst} will  help to carry out this task in the future.

    \noi
Our formulas   are quadratic  in the field strenght and offer
 a starting point for investigating  strong field effects.
In principle, they encompass all kinematic possibilities of the system as a
whole and   permit a description of ultra-relativistic situations, where
$|k \st|^2    \simeq  |k \longi |^2   $.

\noi
In the present state of the art, we notice that, in a {\em weak} field,
 the {\em slow} collective motion (first order in $\epsilon$) of opposite
charges interacting through
a potential which doesnot depend on the energy, escapes the above-mentioned
complication; but in this case the presence of external field results in a
first order Starck effect which obviously vanishes for generic shapes of the
mutual interaction potential.
For the harmonic oscillator for instance, this remark indicates that
 the naive quark model enjoys some kind of stability property.
But if we have perturbation theory in mind,
the computation of significant corrections requires the  setting of a
 nonconventional treatment.

\medskip
\noi
In sofar as approximations are concerned,  it is in order to realize that
two situations are possible:

\noi
Either the magnetic field is considered (like in the previous example) as a
perturbation applied to the system.
Or, in contrast,  the mutual interaction is treated as a perturbation like in
 the helium atom.

\noi
In that latter case, the zeroth order approximation describes two independent
particles  moving in the magnetic field; in this unperturbed motion, the
transverse degrees of freedom are bound by the magnetic field
(corrections to the corresponding spectra are reserved for future work).
We expect to avoid the pathology of "continuous dissolution"
 \cite{dissol}  \cite{brok}  for two reasons:
 The particules we consider here have no spin,
and we  can impose positive individual energies,
requiring that both $P \cdot p_1$ and $P \cdot p_2$ have  positive eigenvalues.

\medskip
\noi The Ansatz which allows for a three-dimensional reduction in our covariant
framework  automatically generates various terms in the wave equation.
We have seen  that the importance of these terms depends on the strenght of
the field and on the  state of motion of the system as  a whole.
Inspection of these terms  indicates that, from a practical point of view,
 the shape of the mutual interaction is "somehow modified"
by the magnetic field.
 As can be read off from  (\ref{4.18}),  the modification implied by
 (\ref{4.16}) is quadratic  in $F$  and may become dominant in  strong fields
provided $V\zer $ is energy-dependent.
This point  concerns most of the realistic two-body potentials.

\noi     Further work is needed  in order to get beyond these qualitative
 indications.

    \noi
For the sake of simplicity we have focused here on scalar particles,
 but naturally  an extension to particles with  spin is desirable.
A generalization to globally charged systems would also be of interest.

\noi
Let us finally mention that, in principle, the contact with more conventional
methods of quantum field theory could be improved, trying to directly derive
 all our terms from  a  Bethe-Salpeter  equation {\em that  takes  the magnetic
 field into account  from the start\/}.       This would
    mean to remake the  work of Bijtebier and Broeckaert \cite{brok} in a
 way which respects the particular symmetry of constant magnetic field, $ i.e.$
    treating all the possible lab frames on the same footing.
To our knowledge, nobody has yet carried out this task.