format: text

text:
# Measurement of Inclusive Spin Structure Functions of the Deuteron

## I Introduction

The nucleon spin structure functions $`g_1^{p,n}(x)`$ and $`g_2^{p,n}(x)`$ and their moments have been extensively studied over the past two decades Ashman et al. (1989); Adams et al. (1997a, b); Anthony et al. (1996); Abe et al. (1998, 1997a); Anthony et al. (2000, 2002); Ackerstaff et al. (1997); Airapetian et al. (1998). At large momentum transfer ($`Q^2>>1`$ (GeV/c)2) and final state mass ($`W>2`$ GeV), these data can be successfully described via perturbative QCD (pQCD) up to next–to–leading order (NLO) and give us access to the helicity-weighted distribution functions $`\mathrm{\Delta }q(x)`$ and $`\mathrm{\Delta }G(x)`$ of quarks and gluons in the nucleon Stratmann and Vogelsang (2001); Gluck et al. (2001); de Florian and Sassot (2000); Abe et al. (1997b). In this kinematic regime, one can relate the first moments $`\mathrm{\Gamma }_1^N=_0^1g_1^N(x)𝑑x`$ of the spin structure functions $`g_1^N(x)`$ ($`N`$ = $`p`$ or $`n`$) to the fraction of the nucleon spin carried by the quark helicities and, via the famous Bjorken sum rule Bjorken (1966, 1970), to the weak axial form factor $`g_A`$.

At lower momentum transfers, $`Q^21`$ (GeV/c)2, corrections proportional to powers of $`1/Q^2`$ develop due to higher twist and target mass effects Stein et al. (1995a, b); Balitsky et al. (1990) in addition to the logarithmic $`Q^2`$ dependence predicted by pQCD. As $`Q^2`$ decreases, an increasing part of the kinematic range $`x=0`$$`1`$ lies in the region of resonant final states ($`W<2`$ GeV), which begin to dominate the spin structure functions. They become less positive (or more negative in the case of the neutron), in particular in the region of the $`\mathrm{\Delta }`$ resonance. Data in this region on structure functions and on the (virtual) photon asymmetries $`A_1`$ and $`A_2`$ for the proton and the neutron,

$`A_1={\displaystyle \frac{\sigma _{1/2}\sigma _{3/2}}{\sigma _{1/2}+\sigma _{3/2}}}={\displaystyle \frac{g_1g_2/\tau }{F_1}},`$
$`A_2={\displaystyle \frac{\sigma _{LT}}{\sigma _{1/2}+\sigma _{3/2}}}={\displaystyle \frac{g_1+g_2}{\sqrt{\tau }F_1}},`$ (1)
can help us unravel the spin-isospin structure of resonance transition amplitudes and their interference with each other and with non-resonant terms. We can also test whether the observed duality between unpolarized deep inelastic and resonant structure functions Bloom and Gilman (1970); Niculescu et al. (2000a, b) is realized for spin structure functions as well Carlson and Mukhopadhyay (1998); Simula et al. (2002). Here, $`\sigma _{1/2}`$ and $`\sigma _{3/2}`$ are the (virtual) photon absorption cross sections for total (photon plus nucleon) helicity 1/2 and 3/2 and $`\sigma _{LT}`$ is the longitudinal–transverse interference cross section, $`F_1`$ is the unpolarized structure function, and $`\tau =\nu ^2/Q^2`$ with $`\nu =EE^{}`$ being the energy loss of the scattered electron.

Due to the dominance of the resonances at low $`Q^2`$, the integrals $`\mathrm{\Gamma }^p`$ and $`\mathrm{\Gamma }^d(\mathrm{\Gamma }^p+\mathrm{\Gamma }^n)/2`$ (which are positive in the scaling region of high $`Q^2`$) decrease rapidly and become negative as $`Q^2`$ approaches zero. In the limit $`Q^20`$, the first moments for the proton and the neutron are constrained by the Gerasimov-Drell-Hearn Sum Rule Gerasimov (1966); Drell and Hearn (1966), which predicts that

$$\mathrm{\Gamma }_1^N(Q^2)\frac{Q^2}{16\pi ^2\alpha }_{\nu _{thr}}^{\mathrm{}}\left(\sigma _{1/2}\sigma _{3/2}\right)\frac{d\nu }{\nu }=\frac{Q^2}{8M^2}\kappa _N^2.$$
(2)
Here, $`\alpha `$ is the fine structure constant and $`M`$ and $`\kappa _N`$ are the mass and anomalous magnetic moment of the nucleon, respectively. Since the GDH sum rule is negative, the integrals $`\mathrm{\Gamma }_1^{p,d}(Q^2)`$ must have a negative slope at $`Q^2=0`$ and then change rapidly at low $`Q^2`$ to meet the positive experimental results in the DIS region.

So far, only phenomenological models for $`\mathrm{\Gamma }_1(Q^2)`$ covering the whole range of $`Q^2`$ exist Anselmino et al. (1989); Soffer and Teryaev (1995, 1997); Burkert and Li (1993); Burkert and Ioffe (1992, 1994); Ji and Unrau (1994). These models are constrained to reach the large–$`Q^2`$ asymptotic value of the integral as measured by deep inelastic data and to approach zero at the photon point with a slope given by the Gerasimov-Drell-Hearn sum rule, Eq. (2). The authors of Refs. Soffer and Teryaev (1995, 1997) use a simple parametrization of the integral $`\mathrm{\Gamma }_{1+2}=[g_1(x)+g_2(x)]𝑑x`$ to interpolate between these two points, and then subtract the integral over $`g_2`$ which is given by the Burkardt-Cottingham sum rule Burkhardt and Cottingham (1970). The approach taken in Refs. Burkert and Li (1993); Burkert and Ioffe (1992, 1994) uses a parametrization of existing resonance data and a vector meson dominance (VMD) inspired interpolation of the remaining integral strength at the two endpoints.

For a complete picture of spin structure functions and their moments, one needs information on both the proton and the neutron. Since free neutron targets are impractical, deuterium (as in the experiment described here) or 3He targets are used instead. An unambiguous extraction of neutron spin structure functions from nuclear ones is less straightforward in the resonance region than in the deep inelastic regime; however, the integrals $`\mathrm{\Gamma }_1^N`$ are much less affected by uncertainties from Fermi motion, off-shell effects, and other nuclear corrections Kaptari et al. (1995); Piller et al. (1996); Ciofi degli Atti et al. (1996a). In particular, studies Ciofi degli Atti et al. (1996b); Arenhovel et al. (1998) show that the integral $`\mathrm{\Gamma }_1^d`$ for the deuteron from pion threshold on up is very close to the incoherent sum of proton and neutron integrals, once a correction for the deuteron D-state has been applied.

So far, only very limited spin structure function data exist in the region of low to moderate $`Q^2`$ and $`W`$ Baum et al. (1980); Abe et al. (1997c), especially on the deuteron. A large program is underway at Jefferson Lab to map out the entire kinematic region $`Q^20.05`$$`5`$ (GeV/c)2 and $`W3`$ GeV. This program consists of measurements on 3He (in Hall A) and on proton and deuteron targets with CLAS (the EG1 collaboration in Hall B). First results from CLAS De Vita et al. (2002) and Hall A Amarian et al. (2002) have already been published.

In the present paper, we present results on the deuteron from the first EG1 run in 1998, in which we measured double spin asymmetries $`A_{||}=D(A_1+\eta A_2)`$ on deuterium with a beam energy of 2.5 GeV. ($`D`$ and $`\eta `$ are kinematical factors, see Sec. IV.) These data cover a range in $`Q^2`$ from 0.27 – 1.3 (GeV/c)2 and final state mass in the resonance region ($`W`$ = 1.08 – 2.0 GeV). The remaining data set from EG1 is presently under analysis and will increase both the kinematic coverage and the statistical precision of our data significantly.

In the following, we give some details on the experiment (Sec. II) and its analysis (Sec. III). We present our results on the deuteron spin asymmetry $`\left(A_1^d+\eta A_2^d\right)(W,Q^2)`$, the structure function $`g_1^d(x,Q^2)`$ and its first moment $`\mathrm{\Gamma }_1^d(Q^2)`$ (Sec. IV), and conclude with a summary and outlook (Sec. V).

## II Experimental Details

The data described in this paper were collected during a three-month run in 1998, as part of the EG1 run group in Jefferson Lab’s Hall B. A polarized electron beam with 2.5 GeV beam energy was scattered off a deuterated ammonia (15ND3) target that was dynamically polarized along the beam direction. The average beam current of 2.5 nA corresponded to an instantaneous luminosity of $`0.4\times 10^{34}`$ cm-2s-1. The beam polarization was measured periodically with a Møller polarimeter and the average beam polarization was 72%.

We used the CEBAF Large Acceptance Spectrometer (CLAS) to detect the scattered electrons. The CLAS detector Mecking et al.  is built around six superconducting coils that produce a toroidal magnetic field. The orientation of the magnetic field can be chosen so that electrons are bent either toward (inbending) or away from the beam line (outbending). The target was placed 55 cm upstream from its normal location in the center of CLAS to lower the angular threshold for electron detection and thus decrease the lower limit on the momentum transfer. Inbending electrons were detected down to a minimum polar angle of 14. During this experiment the geometry of the target excluded particle tracks with a polar angle between 50 and 75. The phi–acceptance is approximately 85%, limited mainly by the torus coils.

The CLAS detector package consists of three layers of drift chambers for track reconstruction, one layer of scintillators for time–of–flight measurements, forward Cerenkov counters for electron–pion discrimination, and electromagnetic calorimeters to identify electrons and neutral particles. A coincidence between the Cerenkov and the calorimeter triggers the data acquisition. Electron particle identification is accomplished using the Cerenkov detector and the distribution of energy deposited in the calorimeter. The large acceptance of CLAS ($``$ 1.5 sr for electrons) and its large kinematic coverage offset the limited luminosity that can typically be reached with polarized solid state targets (of order $`10^{35}`$ cm-2s-1 at best), and allowed us to collect data for the entire $`W`$ and $`Q^2`$ range simultaneously.

The longitudinally polarized target was designed to fit within the 1 m central bore of the CLAS Keith et al. . A pair of superconducting Helmholtz coils provided a 5 T magnetic field along the direction of the electron beam. The magnetic field was uniform to better than $`1\times 10^4`$ in the center of the target over a length of 2 cm and a diameter of 2 cm. The ammonia crystals were contained within a plastic cylindrical cell 1 cm in length and 1.5 cm in diameter. The cell was immersed in a liquid He bath maintained at approximately 1 K by a 4He evaporation refrigerator. The cell was mounted on a target insert that also held a NH3 cell, as well as a 12C and an empty cell. The latter cells were used to study the dilution of the measured asymmetries by events from unpolarized target constituents (see Sec. III.3). The deuterons in the target were polarized using the Dynamic Nuclear Polarization (DNP) technique Abragam and Goldman (1978); Crabb and Meyer (1997) with 140 GHz microwaves. The polarization of the target was monitored on-line using the NMR technique. The NMR results were not used for our final analysis; instead, we extracted the product of beam and target polarization directly from our data, as described in Sec. III.3. The beam was rastered on the ND3 target, although not over the full face of the target. The deuteron polarization suffered from this incomplete raster and from inadequate microwave power and ranged from approximately 10% to 25%. All data were taken with the target polarization along the beam direction, without reversal of the target polarization. The beam helicity was reversed every second.

During the 1998 run, we collected 300 million triggers for an integrated beam charge of about 0.4 mC. From this sample, 100 million electron events passed the cuts described in Sec. III.1. These events covered a kinematic region from the quasielastic region ($`W0.94`$ GeV) to the edge of the deep inelastic region ($`W=2`$ GeV) and for $`Q^2`$ = 0.27 – 1.3 (GeV/c)2. This kinematic coverage is shown in Fig. 1, together with the coverage of the second part of the EG1 experiment.

## III Data analysis

The goal of our analysis was to first determine the double spin asymmetry

$$A_{||}=\frac{\sigma ^{}\sigma ^{}}{\sigma ^{}+\sigma ^{}}$$
(3)
for each kinematic bin and then to extract the physical quantities of interest, the virtual photon asymmetries $`A_1^d+\eta A_2^d`$ and the structure function $`g_1^d`$, from the results. Here, $`\sigma ^{}`$ stands for the differential electron scattering cross section with the target and electron spin pointing in opposite directions along the beam and correspondingly $`\sigma ^{}`$ for parallel target and electron spin.

### III.1 Data selection

For the present analysis, we selected data runs taken with a torus current of +2250 A (inbending electrons) and target polarization parallel to the beam direction. The data were taken with two slightly different beam energies, 2.494 GeV and 2.565 GeV, due to a change of the accelerator configuration. We separated our sample into four different “run groups”, two each with beam energy 2.494 GeV and 2.565 GeV. Each run group corresponds to a contiguous set of runs with the same target material and approximately constant target parameters and running conditions. Only runs with stable beam and detector performance were included in our sample. The 2.565 GeV groups also contained carbon target runs that were used to determine the dilution factor (see Sec. III.3). We analyzed events with scattering angles from about 14 to 50 and scattered electron energies from 0.5 GeV to 2.5 GeV.

The data were sorted according to the helicity of the electron beam. During our run, the beam helicity followed a “pseudo–random” pattern of helicity pairs, where the first “bucket” (of 1 second length) of each pair was given random helicity and the second its complement. We matched the sequence of helicity bits for each event with the pattern sequence recorded in helicity scalers and discarded pairs for which the helicity assignment was inconsistent. We also discarded pairs with significantly different (by more than 10%) beam intensity in the two buckets (due to beam fluctuations or trips). The final data sample contained only matched pairs of buckets with stable running conditions.

All events were accumulated in small bins of $`W`$ ($`\mathrm{\Delta }W=0.02`$ GeV) and $`Q^2`$ ($`\mathrm{\Delta }Q^2/Q^220\%`$), separately for both beam helicities. (The data on asymmetries and $`g_1`$, shown in Sec. IV, are weighted averages of several such bins). In addition, we also accumulated the integrated beam charge for each of the helicity buckets (corrected for deadtime) to normalize the helicity-sorted counts in each bin. We found that on average, there was a 0.3% difference between the integrated charge for the two opposite helicities, possibly stemming from the sensitivity of the photocathode in the polarized source to small remaining linear polarization components or beam motions of the photo–ionization laser beam. Our normalization method removed the effect of this asymmetry, and it was further suppressed by reversing the relative sign between the helicity at the cathode and in the experimental Hall (through spin–precession in the injector and the accelerator).

### III.2 Electron cuts

We selected electron events by first requiring a negative track with matching signals in the time–of–flight (ToF) scintillators, the Cerenkov counters (CC), and the electromagnetic calorimeter (EC). In the presence of several such tracks, the track with the shortest flight time was selected as the electron candidate. Some additional cuts on the track vertex along the beam line removed events from the entrance and exit windows of the polarized target chamber, as well as badly reconstructed tracks.

We used information from the CC and the EC to further separate electrons from negative pions. We required a signal in the CC that exceeded 50% of the average signal for a single photo–electron. Furthermore, we required that the energy measured in the EC exceeded 20% of the candidate electron momentum (the average sampling fraction of the EC was 27%). A typical example for the ratio of sampled EC energy over momentum is shown in Fig. 2. The open histogram shows events that passed all other electron cuts (including the CC cut).

We also collected a sample of $`\pi ^{}`$ events with no signal above threshold in the CC. As shown by the shaded area in Fig. 2, the $`E/p`$ spectrum associated with $`\pi ^{}`$ events is strikingly different from the electron spectrum. Under the conservative assumption that all events below a $`E/p`$ ratio of 0.15 came from pions, we cross–normalized the two spectra below that point and estimated the remaining pion contamination of our electron sample by the ratio of the two integrated spectra above our cut of $`E/p>0.2`$. For all kinematics studied, this remaining contamination turned out to be less than 1%.

The reconstructed momenta of the scattered electrons were corrected for effects from unknown torus field distortions and slight drift chamber misalignments. We used NH3 runs taken interleaved with the ND3 ones to determine the correction factor by optimizing the position and width of the elastic peak ($`W=0.938`$ GeV) for all scattering angles $`\theta `$ and $`\varphi `$. The resulting corrections were of the order 0.1% on average.

### III.3 Dilution and polarization

The double spin asymmetry $`A_{||}`$ can be extracted from the count rate asymmetry (normalized by the integrated beam charge) after accounting for the dilution from unpolarized target constituents and the beam ($`P_b`$) and target ($`P_t`$) polarization:

$$A_{||}^{meas}=\frac{1}{DFP_bP_t}\frac{N^+/Q^+N^{}/Q^{}}{N^+/Q^++N^{}/Q^{}},$$
(4)
where $`N^{+,}`$ are the counts and $`Q^{+,}`$ are the integrated beam charge for positive and negative helicity.

We determined the dilution factor $`DF`$ in Eq. (4) by approximating the contribution to the count rates from unpolarized target constituents (target foils, LHe coolant and 15N in ammonia) with the spectra taken on the carbon target. Some components of these two targets were the same (e.g., the LHe coolant and foils were present for the carbon target as well), and carbon, nitrogen, and even 4He have similar binding energies per nucleon and Fermi momenta, suggesting that their inclusive electron scattering spectra are similar after correcting for the total number of target nucleons. (This assumption has since been verified to better than 3% with dedicated runs on a pure 15N target during the second part of EG1).

To account for the different number of nucleons in each target and different overall target thicknesses, we cross–normalized the carbon target spectra to the ammonia target spectra. We determined a normalization constant $`A`$ such that the two spectra had the same number of counts below a cut–off missing mass $`W_{cut}`$, well below the quasi–elastic peak. The cut–off ranged from $`W_{cut}=0.835`$ GeV at $`Q^2=0.3`$ (GeV/c)2 to $`W_{cut}=0.5`$ GeV at $`Q^2=1.2`$ (GeV/c)2, and was chosen so that the deuteron contribution was negligible, according to a Monte Carlo simulation of the deuteron wave function.

The dilution factor can then be written as

$$DF=\frac{N_{ND_3}AN_C}{N_{ND_3}},$$
(5)
where the numerator is the count rate due to deuterium alone. The results of this method for an intermediate $`Q^2`$ bin are shown in Fig. 3. The normalized carbon spectrum (circles) has been subtracted from the ammonia spectrum (solid triangles) to yield the deuteron spectrum (open triangles). The line indicates the result of our Monte Carlo simulation of the deuteron spectrum alone, which is based on quasi-elastic scattering (PWIA) and the Paris wave function Lacombe et al. (1980) for the deuteron. The dilution factor for our experiment was around $`DF0.2`$.

The second ingredient needed in Eq. (4) is the product of beam and target polarization. We measured both the beam polarization (with a Møller polarimeter) and the target polarization (using NMR) individually during the run. However, due to the small amount of target material and its inhomogeneous exposure to the electron beam, the NMR results were not very precise and reliable. Instead, we determined directly the product $`P_bP_t`$ by extracting it from the measured asymmetry in the quasi–elastic region. For this purpose, we used inclusive quasi–elastic events $`d(e,e^{})`$ in the range 0.85 GeV$`W`$ 1.0 GeV.

The asymmetry $`A_{||}`$ for elastic scattering from protons and neutrons can be calculated from known nucleon form factors with very little systematic uncertainty (less than 1-2% in our kinematic region). We used our simulation of the deuteron wave function to calculate the expected asymmetry for inclusive quasi–elastic scattering within our kinematic cuts, which differed only slightly from the cross section–weighted average of the proton and neutron asymmetries. We used the dilution factor determined by the method described above to extract the product $`P_bP_t`$.

Due to the large kinematic coverage of CLAS, data on the quasi–elastic asymmetries were collected continuously and simultaneously with the inelastic asymmetry data. The extracted average polarization product $`P_bP_t`$ for each of the four run groups is therefore a faithful representation of the running conditions for that group, with minimal systematic uncertainties. Our results are shown in Fig. 4, where we divided the product $`P_bP_t`$ by the measured beam polarization to extract the target polarization. The results for each of the individual run groups have statistical errors on the order of 13%, which were included in the total statistical error of the asymmetries from each run group. The final results for the inelastic asymmetries are statistically weighted averages from the four run groups, with a contribution to their statistical errors from the polarization product of about 6.7% of their values.

### III.4 Other backgrounds

After dividing out the dilution factor and the beam and target polarizations in Eq. (4), we corrected the extracted asymmetry for additional background contributions. These include contamination of the scattered electron sample by negative pions and pair-produced electrons, as well as contributions from polarized target constituents other than deuterium.

We already discussed the contribution from pions misidentified as electrons, which was less than 1% in all cases. A more important contribution comes from electrons that are decay products of neutral pions (either through the Dalitz decay $`\pi ^0\gamma e^+e^{}`$ or pair conversion of decay photons). The rate of electrons from these decays was estimated using the Wiser fit Wiser (1977) for pion photoproduction and tested against the Monte Carlo code “PYTHIA”. We also measured directly the rate of positron production in each kinematic bin (again making use of the large acceptance of CLAS for both positively and negatively charged particles). This rate should be equal to that of electrons from charge-symmetric decays and was found to agree well with the Wiser fit. The asymmetry for positrons was found to be consistent with zero and in any case no larger than the asymmetry for electron scattering events. We used a parametrization of our results to estimate the fraction of detected electrons coming from these decays. This fraction was typically 1% for most of the kinematic region, but increased up to 20% at the highest $`W`$ values. We corrected our data for this background by applying a further dilution factor to our asymmetries. Since we could not exclude a small non–zero asymmetry for these events, we assumed a systematic uncertainty equal to the size of this correction.

The nitrogen in our dynamically polarized ammonia target carries a small residual polarization, which leads to a partially polarized bound proton in 15N. Possible additional polarized target species include isotopic impurities of 14N and 1H. Extensive experience with similar targets at SLAC Abe et al. (1998) shows that the corresponding corrections to the asymmetry are at most a few percent. We included the uncertainty due to these contributions in our systematic error.

Another potential contribution to the measured asymmetry comes from parity-violating electron scattering off all target constituents. However, at the low momentum transfers of our experiment, the expected asymmetry is less than $`10^4`$ Prescott et al. (1979) and can be treated as another (small) systematic uncertainty.

### III.5 Radiative corrections and models

The final step in the extraction of the desired “Born asymmetry” $`A_{||}`$ requires correcting the measured asymmetry for higher–order electromagnetic processes (internal radiative corrections) and electron energy loss through bremsstrahlung in the target before or after the scattering (external radiative corrections). These radiative corrections were applied separately to the numerator and the denominator of Eq. (3), which yields an additive ($`A_{RC}`$) and a multiplicative ($`F_{RC}`$) correction term:

$$A_{||}=A_{||}^{meas}/F_{RC}+A_{RC}.$$
(6)
Here, the factor $`1/F_{RC}`$ represents the increase of the denominator in Eq. (3) due to the radiative elastic and quasi-elastic tails that act like an additional dilution of the inelastic events. Correspondingly, the statistical error of the final result was scaled up by $`1/F_{RC}`$ as well.

Both components ($`F_{RC}`$, $`A_{RC}`$) were determined by running the code “RCSLACPOL” developed at SLAC Abe et al. (1998). This code uses parametrizations of all relevant input quantities (structure functions and form factors), as well as a model of our target, to calculate both fully radiated and Born cross sections and asymmetries. It is based on the approach developed by Kuchto and Shumeiko Kukhto and Shumeiko (1983) for the internal corrections and by Tsai Tsai (1974) for the external corrections, including the radiative depolarization of the beam due to external bremsstrahlung.

We used parametrizations of the world data on polarized and unpolarized structure functions and elastic form factors as input for the radiative correction code and to extract physics quantities of interest from the measured asymmetries. These parametrizations are described in Anthony et al. (2000) and are based on fits to unpolarized structure function data from NMC Arneodo et al. (1995) and SLAC Bodek et al. (1979); Tao et al. (1996); Stuart et al. (1998); Abe et al. (1999) and polarized structure function data from SLAC Baum et al. (1980); Anthony et al. (1996); Abe et al. (1998, 1997a); Anthony et al. (2000, 2002), CERN Ashman et al. (1989); Adams et al. (1997a, b), and HERMES Ackerstaff et al. (1997); Airapetian et al. (1998). The nucleon form factors were taken from Ref. Bosted (1995) with updated values for the ratio $`G_{Ep}/G_{Mp}`$ from the recent Jefferson Lab experiment Jones et al. (2000). For the asymmetries $`A_1`$ and $`A_2`$ in the resonance region, we used parametrizations of resonance transition amplitudes from Ref. Burkert and Li (1993) (in the form of a computer code named “AO”) and Ref. Drechsel et al. (1999) (“MAID”) together with a fit of the SLAC data Abe et al. (1997c). We also included our own preliminary asymmetry data in these fits. All fits were varied within reasonable errors or replaced with alternative existing fits to study the systematic dependence of our final results on these parametrizations.

### III.6 Systematic errors

The total systematic error on our data ranges from 25% to 50% of the statistical error for the asymmetries and from 35% to 50% of the statistical error for the structure function $`g_1^d`$. The leading contributions to these systematic uncertainties come from radiative corrections (40–50% of the total systematic error on average), uncertainties in the unpolarized structure functions needed to extract final physics results (also 40–50% of the total), and the dilution factor (about 40%). We also considered the effect of finite resolution and errors in the measured kinematic variables (about 10% of the total). At higher $`Q^2`$ and especially higher $`W`$, pair-symmetric decay electrons also contributed significantly to the overall systematic uncertainty (15–20% averaged over all kinematic bins and most of the systematic error at the kinematic limit). Finally, for the extraction of the spin structure function $`g_1^d`$ and its integrals, some model assumption about the virtual photon asymmetry $`A_2`$ is needed (see Sec. IV) and leads to a further systematic error (up to 50%).

We accounted for each of these systematic errors by changing a relevant input parameter or model, and then repeating the entire analysis up to the final results, including the integrals of $`g_1^d`$ over the measured region. We took the error as the deviation of the alternative results from the standard analysis. We added all uncorrelated systematic errors in quadrature. The final systematic errors are shown in the data tables in Sec. IV.

For the radiative correction errors, we varied all input models and parametrizations for the radiative code, including polarized and unpolarized structure functions, form factors, and the target model, within realistic limits. We also checked the accuracy of the peaking approximation by comparing the results with those from a full integration without approximations.

Similarly, we varied the models for the unpolarized structure functions $`F_1^d`$ and $`R=\frac{\sigma _L}{\sigma _T}`$, which entered the extraction of $`g_1^d`$ and the asymmetry $`A_1^d+\eta A_2^d`$ from our data (see Sec. IV). We used different fits of the world data Bodek et al. (1979); Abe et al. (1999); Ricco et al. (1999), and studied their effect on the final physics results. In the case of the polarized structure function $`g_1^d`$ and its integrals, we also varied the model for the asymmetry $`A_2^d`$ from $`A_2^d=0`$ to the prediction by the “MAID” code and a simple parameterization based on the twist–2 result by Wandzura and Wilczek Wandzura and Wilczek (1977) that describes the SLAC data Anthony et al. (2002) well.

For the error introduced by the uncertainty in the dilution factor, we varied the cross–normalization between the carbon and ammonia target data by an amount of 6%, consistent with the variations observed for different $`W`$ and $`Q^2`$ ranges and possible differences in the 12C and 15N spectra. This yields an average variation of the dilution factor Eq. (5) of 25%, making this error a safe upper bound for all systematic errors that are directly proportional to the measured asymmetry.

The CLAS momentum resolution and reconstruction effects were studied by moving all data points by 0.02 GeV in $`W`$ and recalculating the final results. The effect of this variation on the integrals of $`g_1^d`$ also gave an upper limit to systematic errors due to the integration method, which consisted in a simple sum of all bins of size $`\mathrm{\Delta }W=0.02`$ GeV, multiplied by the bin width in $`x`$.

Other systematic errors were either negligible or have already been described in the previous section. We note that we do not have a significant systematic error from the beam and target polarization product, since they were directly determined from our data (with minimal theoretical uncertainty). In particular, the theoretical asymmetry $`A_{||}^{elas}`$ is only weakly dependent ($`\pm 1\%`$) on the elastic form factor ratio $`G_E/G_M`$ for the proton. However, the statistical error of this method is not negligible and was included in the total statistical error of the final results.

## IV Results

### IV.1 Virtual photon asymmetries

We extracted a combination of the virtual photon asymmetries, $`A_1^d+\eta A_2^d`$, from our data on $`A_{||}`$ using a parametrization Abe et al. (1999) of the structure function $`R`$, via the relationship

$$A_{||}=D\left(A_1+\eta A_2\right),$$
(7)
where the virtual photon depolarization factor is given by $`D=(1ϵE^{}/E)/(1+ϵR)`$ and $`\eta =ϵ\sqrt{Q^2}/(EϵE^{})`$ ($`ϵ`$ is the virtual photon polarization parameter, $`E`$ is the beam energy and $`E^{}`$ is the scattered electron energy).

The extracted photon asymmetries $`(A_1^d+\eta A_2^d)(W,Q^2)`$ for three different $`Q^2`$ bins are listed in Tables 13, together with their statistical and full systematic errors. We show the results for our intermediate $`Q^2`$ bin in Fig. 5, together with previous data from SLAC Abe et al. (1998) and some model calculations. A comparison of the three different $`Q^2`$ bins can be found in Fig. 6.

Since we didn’t measure the asymmetry with the target polarization perpendicular to the electron beam ($`A_{}`$), we cannot directly extract the asymmetry $`A_1^d`$ or $`A_2^d`$. The interference term $`A_2`$ is limited by $`|A_2|<\sqrt{R(A_1+1)/2}`$, where the value of $`R`$ is around 0.1 - 0.3 at $`Q^2=0.5`$ (GeV/c)2 Abe et al. (1999) and the typical size of $`\eta `$ for our experiment ranges from 0.1 at $`W=2`$ GeV to 1.2 right at pion threshold ($`W=1.08`$ GeV). Correspondingly, the asymmetry $`A_2^d`$ could contribute as much as 0.07 (high $`W`$) to 0.15 (at threshold) to the asymmetries shown in Figs. 5 and 6. However, according to our parametrization, this contribution should be more typically of order 0.02.

With this caveat, one can conclude that the data shown in Fig. 5 exhibit the expected behavior for the asymmetry $`A_1^d`$. In the region of the $`\mathrm{\Delta }`$(1232) resonance, the asymmetry is strongly negative, and fully compatible with the expectation $`A_1^d=0.5`$ for the resonance contribution alone. Beyond $`W=1.4`$ GeV, the asymmetry becomes positive, indicating that helicity–1/2 transition amplitudes begin to dominate even at this rather low $`Q^2`$. However, even in the region of the S11 resonance the asymmetry is markedly smaller (around 0.15) than for the proton (around 0.5, see Ref. Abe et al. (1998)), indicating that for the neutron alone the helicity–3/2 amplitude may still be larger. Figure 5 also shows the predicted full asymmetry from our parametrization and a prediction for the resonance contributions to $`A_1^d`$ alone. The latter is based on the code “AO” Burkert and Li (1993), which uses a fit of exclusive pion electro- and photoproduction data to parametrize resonant and Born pion production amplitudes. Apparently, the contribution from the resonances alone already describes the data well in the region of low to intermediate $`W`$, while non–resonant contributions (and maybe a sizable asymmetry $`A_2^d`$) are needed at high $`W`$. In general, our data agree fairly well with model predictions and the existing SLAC data. However, they have significantly smaller statistical errors and better resolution in $`W`$, as well as coverage down to lower $`Q^2`$ than the SLAC data.

A comparison of our results for different $`Q^2`$ (see Fig. 6) shows a general trend toward more positive asymmetries for higher $`Q^2`$, especially in the region of the S11 and D11 resonances. This is in agreement with the expected transition from helicity–3/2 dominance at low $`Q^2`$ (and especially at the photon point, where it yields the negative value for the GDH sum rule), and helicity–1/2 dominance at higher $`Q^2`$. In the limit of very large $`Q^2`$, the asymmetry $`A_1^d`$ in the resonance region should become close to 1, as predicted by pQCD as well as hyperfine-improved quark models and duality arguments. A similar behavior is observed for the proton asymmetries Abe et al. (1998).

### IV.2 Spin structure function $`g_1^d`$

The spin structure function $`g_1^d(W,Q^2)`$ was calculated from the photon asymmetry $`(A_1^d+\eta A_2^d)(W,Q^2)`$ for each bin using

$$\begin{array}{c}g_1^d(W,Q^2)=\frac{\tau }{1+\tau }\left(A_1^d+\frac{1}{\sqrt{\tau }}A_2^d\right)F_1^d(W,Q^2)\hfill \\ =\frac{\tau }{1+\tau }\left(\left(A_1^d+\eta A_2^d\right)+\left(\frac{1}{\sqrt{\tau }}\eta \right)A_2^d\right)F_1^d(W,Q^2).\hfill \end{array}$$
(8)
Here, $`F_1^d(F_1^p+F_1^n)/2`$ represents the unpolarized structure function of the deuteron (per nucleon) and $`\tau =\nu ^2/Q^2`$. Because of the partial cancellation of the two terms in $`\left(\frac{1}{\sqrt{\tau }}\eta \right)`$, $`g_1^d`$ is less sensitive to the asymmetry $`A_2`$. We list our results for $`g_1^d`$ with their statistical and full systematic errors (including the uncertainty due to $`A_2`$) in Tables 4 \- 6.

In Fig. 7, we show our results for all three values of $`Q^2`$, plotted against the Nachtmann scaling variable $`\xi =Q^2/M(\nu +q)`$. This variable corresponds to Bjorken $`x`$ at high $`Q^2`$ while it takes target nucleon mass corrections into account and therefore reduces “kinematical higher twist” scaling-violating effects at lower $`Q^2`$. Together with our data, we also show as reference the prediction for $`g_1^d(\xi ,Q^2=5`$ (GeV/c)$`{}_{}{}^{2})`$ from our model. The assumption of local quark-hadron duality predicts that structure functions like $`F_1`$ and $`g_1`$ should, on average, approach a universal scaling curve if plotted versus the variable $`\xi `$, even in the resonance region. This is confirmed down to rather low $`Q^2`$ in the case of the unpolarized structure function $`F_2^p`$ Niculescu et al. (2000a, b). Apparently, local duality does not work as well for the polarized structure function $`g_1^d`$ at high values of $`\xi `$ where the asymmetry is dominated by the $`\mathrm{\Delta }`$ resonance and therefore is negative. Overall, the approach to the “asymptotic value” for $`Q^2=5`$ (GeV/c)2 seems to be relatively slow; only our highest $`Q^2`$ bin shows fairly good agreement beyond the region of the $`\mathrm{\Delta }`$ resonance.

### IV.3 Integrals

We calculated the integrals $`\mathrm{\Gamma }_1^d(Q^2)=g_1^d(x,Q^2)𝑑x`$ for our results on $`g_1^d(x,Q^2)`$ over the (ordinary) Bjorken variable $`x`$ for four different $`Q^2`$ bins, beginning at quasi–free pion production threshold ($`W`$ = 1.08 GeV) up to the kinematic limit of our data. (The first two $`Q^2`$ bins are the same as shown in Tables 4 and 5, while we split the last bin into two halves, from $`Q^2`$ = 0.65 to 0.92 (GeV/c)2 and from $`Q^2`$ = 0.92 to 1.3 (GeV/c)2.) We expect that these integrals are close to an incoherent average over the individual nucleons (proton and neutron) in deuterium, reduced by the D-state correction factor $`(11.5P_D)`$, where $`P_D0.05`$ is the deuteron D-state probability. The results are shown in the third column of Table 7 and the upper kinematic limits for $`W`$ are listed in the second column. These upper $`W`$ bounds correspond to lower limits of $`x=(0.1,0.15,0.21,0.32)`$ for the four $`Q^2`$ bins, respectively.

We use our model to estimate the contribution to the integral below these limits and show the resulting “full” integrals and their systematic errors in the last two columns of Table 7. These systematic errors include a contribution from the uncertainty of this extrapolation to $`x=0`$. To estimate this uncertainty, we studied the variation of the low–$`x`$ contribution according to different fits to the world data; also, since there are few high–precision data below $`x=0.03`$, we added a systematic error equal to the value of the integral from $`x=0`$ up to 0.03. Due to the large theoretical uncertainty about the shape of the spin structure functions at very low $`x`$ and the absence of high-precision data in this region, the error on this extrapolation may be even larger than indicated by our systematic error estimate (see below).

Our results for the first moment $`\mathrm{\Gamma }_1^d(Q^2)`$ of the spin structure function $`g_1^d`$ are shown in Fig. 8. The solid line at higher $`Q^2`$ is a fit to the world’s data in the DIS region including QCD corrections up to second power in the strong coupling constant. The dotted line indicates the slope for the integral at $`Q^2=0`$ predicted by the GDH sum rule (we use the incoherent sum of the results for the proton and for the neutron, normalized to two). The short–dashed line is the result from the code “AO” Burkert and Li (1993) for the contribution from the nucleon resonances only. The long–dashed line by Burkert and Ioffe Burkert and Ioffe (1992, 1994) is the “AO” result plus a term that depends smoothly on $`Q^2`$ and interpolates between the part that is missing at $`Q^2=0`$ to saturate the GDH sum rule and the full value of $`\mathrm{\Gamma }_1`$ in the high $`Q^2`$ limit. Fig. 8 also shows the prediction from the model by Soffer and Teryaev Soffer and Teryaev (1995, 1997) (dot–dashed line). They use an interpolation of the integral over the structure function $`g_T=g_1+g_2`$, which converges to $`\mathrm{\Gamma }_1`$ at high $`Q^2`$ and remains positive down to the photon point where its slope is given by a combination of the nucleon charge and anomalous magnetic moment. They subtract the contribution from the integral over $`g_2`$ (which is related to nucleon form factors via the Burkardt–Cottingham sum rule) to obtain the integral $`\mathrm{\Gamma }_1`$ alone. The same authors have recently published a new parametrization of the proton–neutron difference integral for all $`Q^2`$ Soffer and Teryaev (2002) which might change the curve for the deuteron shown here. The solid triangles are based on EG1 data alone and the open triangles include the estimated contribution to the integral from beyond our kinematic limits. The inner error bars are statistical and the outer error bars represent the systematic errors added in quadrature. They include the uncertainty on the estimated low-$`x`$ contribution for the full integrals (open triangles).

The first conclusion one can draw from Fig. 8 is that the integral over our measured region (essentially the resonance region) is in rather good agreement with the prediction of the “AO” parametrization for resonance contributions only. The data follow the predicted trend from negative values at small $`Q^2`$, where the $`\mathrm{\Delta }`$ resonance contributes most of the integral and most other resonances are also dominated by the helicity–3/2 transition amplitude, to positive values at higher $`Q^2`$, where the helicity–1/2 amplitude begins to take over and the importance of the $`\mathrm{\Delta }`$ is diminished. Since we did not include Born terms or other non-resonant terms in the curve labeled “AO”, one can conclude that these terms must contribute relatively little to the integral over the resonance region in the case of the deuteron. This may be due to a partial cancellation between the asymmetry of the proton (which is likely positive for these terms) and that of the neutron.

Extrapolating the integral down to $`x=0`$ seems to change the results only moderately (in the negative direction at low $`Q^2`$ and towards more positive values at higher $`Q^2`$). This can be understood again as a cancellation between a strongly negative-going trend of the structure function $`g_1^n(x)`$ as $`x`$ goes to zero and a more positive trend for $`g_1^p(x)`$, according to existing DIS data and NLO analyses Abe et al. (1997b); Anthony et al. (1999). However, at present our understanding of the behavior of spin structure functions at very low $`x`$ is still incomplete, making this extrapolation rather uncertain (as it is in the DIS region). Therefore, the error bars on our open triangles may still underestimate that uncertainty. The emergence of new information on the low–$`x`$ behavior of spin structure functions over the past five years is responsible for most of the apparent disagreement between our quoted results and those from the E143 experiment at SLAC. The integrals over the resonance region alone agree fairly well with the SLAC data (to within 1.1 standard deviations); however, the extrapolation beyond $`W=2`$ GeV is much more negative for the parametrization used in the present analysis and would move the SLAC data points down by about 0.008 and 0.015 at $`Q^2=0.5`$ and 1.1 (GeV/c)2, respectively. With this proviso, our data are (marginally) consistent with the SLAC data, but have much improved statistical errors and cover lower $`Q^2`$.

Our data lie somewhat below both phenomenological predictions for the full integral shown in Fig. 8, suggesting a slower transition from the negative values near the photon point to the positive asymptotic value at high $`Q^2`$. The zero–crossing appears to occur somewhere between $`Q^2=0.5`$ (GeV/c)2 and $`Q^2=0.8`$ (GeV/c)2, significantly later than in the case of the proton Abe et al. (1998). However, the systematic errors are highly correlated point–to–point so that the deviation from the predictions by Burkert and Ioffe Burkert and Ioffe (1992, 1994) and by Soffer and Teryaev Soffer and Teryaev (1995, 1997) is not highly significant.

## V Summary and Outlook

In this paper, we report first results on deuterium for inclusive spin structure functions in the nucleon resonance region from the EG1 program at Jefferson Lab. These data significantly expand the kinematic coverage and statistical precision beyond the only previous data from SLAC Abe et al. (1997c). We find generally reasonable agreement between these two data sets and various model predictions and parametrizations. In particular, the importance and the negative asymmetry of the $`\mathrm{\Delta }`$ resonance is confirmed, as is the general trend to more positive asymmetries at higher $`Q^2`$ and $`W`$.

The spin structure function $`g_1^d`$ is less positive than for the proton case Abe et al. (1998), indicating that the neutron contribution is mostly negative in our kinematic region. While $`g_1^d(\xi ,Q^2)`$ seems to be approaching the DIS scaling curve for large $`W`$ and $`Q^2`$, there are significant deviations from “local duality”, again mostly due to the $`\mathrm{\Delta }`$ resonance.

The integral over $`g_1^d`$ follows the expected trend in general, rising towards the DIS limit at the highest measured $`Q^2`$ while dropping rapidly below 0 towards our lowest $`Q^2`$ point. Clearly, neither the kinematic reach (in $`W`$ and $`Q^2`$) nor the statistical precision of the present data set allow a definite statement about the validity of (or the approach towards) the GDH sum rule limit. However, our data constrain the general trend required of any theory that aims to describe the spin structure of the nucleon over the full range of length scales, from the real photon point to the scaling limit.

Spin structure function data on the deuteron, together with the corresponding proton results, should in principle allow us to separate the different isospin contributions to the resonant and non-resonant asymmetries. However, the first run of EG1 analyzed here did not yield enough statistical precision to make a direct separation of proton and neutron contributions to the deuteron asymmetry feasible. However, we plan to submit results on the integral $`\mathrm{\Gamma }_1`$ for the neutron and the proton–neutron difference, extracted from our data on the proton and the deuteron, in a separate paper. In the meantime, the complete EG1 data set has been collected in a second run, which will yield a nearly tenfold improvement in statistics for the deuteron and a wider coverage towards both lower and higher $`Q^2`$ and higher $`W`$. Once analyzed, this vastly larger data set will allow us to investigate in detail resonance electro–production on the neutron and the approach of the first moment of $`g_1^d`$ and $`g_1^n`$ towards the GDH sum rule at the real photon point.

###### Acknowledgements.

We would like to acknowledge the outstanding effort of the Accelerator, Target Group, and Physics Division staff at TJNAF that made this experiment possible. This work was supported by the U.S. Department of Energy, the Italian Istituto Nazionale di Fisica Nucleare, the U.S. National Science Foundation, the French Commissariat à l’Energie Atomique, the French Centre National de la Recherche Scientifique, and the Korea Science and Engineering Foundation. The Southeastern Universities Research Association (SURA) operates the Thomas Jefferson National Accelerator Facility for the United States Department of Energy under contract DE-AC05-84ER40150.


text:
# Non-ideality of quantum operations with the electron spin of a 31P donor in a Si crystal due to interaction with a nuclear spin system.

also at \]Physics Department, Kazan State University, Kazan 420008, Russia

## Abstract

We examine a 31P donor electron spin in a Si crystal to be used for the purposes of quantum computation. The interaction with an uncontrolled system of 29Si nuclear spins influences the electron spin dynamics appreciably. The hyperfine field at the 29Si nuclei positions is non-collinear with the external magnetic field. Quantum operations with the electron wave function, i.e. using magnetic field pulses or electrical gates, change the orientation of hyperfine field and disturb the nuclear spin system. This disturbance produces a deviation of the electron spin qubit from an ideal state, at a short time scale in comparison with the nuclear spin diffusion time. For $`\mathrm{H}_{\mathrm{ext}}9`$ T, the estimated error rate is comparable to the threshold value required by the quantum error correction algorithms. The rate is lower at higher external magnetic fields.

The interest in the evolution of a single electron spin interacting with a system of uncontrolled nuclear spins (or spin bath) arose from proposed models of spin qubits and qubit-qubit manipulations Loss ; Kventsel ; Kane ; Awschalom ; Jiang ; Hu ; Band ; Glasser . The electrons confined in a quantum dot or in a hydrogen state around a shallow donor impurity are influenced by many nuclear spins through the hyperfine dipole-dipole interaction. At low temperatures, where spin-phonon decoherence Mozyrsky  is ineffective, the dominant mechanism for the electron spin decoherence is via this interaction Khaetskii ; Efros ; Saykin ; Sousa . Additional error in a qubit state appears from non-ideal quantum gate operations XHu . In this work we consider the effect of the nuclear spin bath on results of the quantum operations.

The approximation of short bath correlation times and weak system-bath interaction, which leads to a Markovian-like equation for relaxation and decoherence, is hardly applicable for quantum computation. In that case the evolution of qubit at short times is important Privman . Different terms of system-bath interaction and the bath dynamics influence the evolution of a single qubit at different time scales. In particular, for the spin qubit, interacting with a spin bath, several time scales can be defined Stamp .

The irreversible dynamics of the nuclear spin bath determines the longest time scale. At low temperatures this can be explained in terms of spin diffusion Slichter . The characteristic diffusion time is determined by the flip-flop transitions ($`I_+^1I_{}^2+I_{}^1I_+^2`$) for two neighbouring nuclear spins. These processes are much slower in comparison with the electron spin dynamics. For shorter time intervals, the electron spin sees a ”frozen” system of nuclear spins. The diagonal component, $`S_zI_z`$, of the hyperfine interaction with a ”frozen” spin system leads to the constant shift of the electron Larmor frequency. It is possible to eliminate this in quantum algorithms by various additional correction operations (in analogy with a $`\pi /2\pi `$ pulse sequence), or by the preliminary measuring of the ”nuclear magnetic field” at the localized electron. The single transitions with the electron spin flip ($`S_+I_\alpha +S_{}I_\alpha `$, $`\alpha =+,,z`$) are suppressed by energy conservation, unless the external magnetic field is sufficiently low Khaetskii ; Saykin . However for nuclear spin flips, this suppression is not as strong, due to the small Zeeman splitting of nuclear levels.

The effective magnetic field at each nuclear position is formed by the external and hyperfine fields. It is non-colinear to the external field, due to the hyperfine terms, $`S_zI_++S_zI_{}`$. Quantum operations with the electron spin reorient this field, and cause the nuclear spin bath to oscillate near the nuclear Larmor frequency. These oscillations can be appreciably faster compared to nuclear spin diffusion. As a result, this produces a deviation of the electron spin qubit from an ideal state. We study this effect in detail in the case of a 31P donor electron spin qubit system.

If, initially, at the time $`t=0^{}`$ the qubit is uncoupled with the nuclear spin bath, and the bath is in an equilibrium state, $`\rho _\sigma (0^{})\rho _{bath}(0^{})`$, after very fast quantum operation the system-bath state will be $`\rho _\sigma (0^+)\rho _{bath}(0^{})`$. We assume here, that the qubit operation does not change the bath state significantly. The qubit density matrix at the time $`t`$ is

$$\rho _\sigma (t)=\mathrm{Tr}_{bath}\left\{e^{\frac{i}{\mathrm{}}Ht}\rho _\sigma (0^+)\rho _{bath}(0^{})e^{\frac{i}{\mathrm{}}Ht}\right\}.$$
(1)
At short times, in comparison with the nuclear spin diffusion, we neglect by the nuclear dipole-dipole interaction. The magnetic Hamiltonian of a single electron spin localized at a phosphorus donor in a 29Six 28Si1-x crystal will be Hale

$$H=g\mu _\mathrm{B}\mathrm{𝐇𝐒}+a_0\mathrm{𝐒𝐈}_0\gamma _\mathrm{P}\mathrm{}\mathrm{𝐇𝐈}_0+\underset{k}{}\left(\mathrm{𝐒𝐀}(𝐑_k)𝐈_k\gamma _{\mathrm{Si}}\mathrm{}\mathrm{𝐇𝐈}_k\right),$$
(2)
where $`𝐈_0`$ is a donor nuclear spin and $`𝐈_k`$ are 29Si nuclear spins. The electron spin interaction with the phosphorus nucleus, the second term on the right-hand side of Eq. (2), is isotropic. The tensor, $`𝐀(𝐑_k)`$, describes the hyperfine interaction of an electron spin with a nuclear spin of 29Si at a lattice position $`𝐑_k`$. It consists of isotropic and anisotropic parts as

$$A_{\alpha \beta }(𝐑_k)=\frac{8}{3}\pi g_0\gamma _{\mathrm{Si}}\mathrm{}\mu _\mathrm{B}|\psi (𝐑_k)|^2\delta _{\alpha \beta }+D_{\alpha \beta }(𝐑_k),$$
(3)
where we use the free electron $`g`$-factor Paget . Values of hyperfine constants have been measured by G. Feher Feher  and later by E. B. Hale and R. L. Mieher Hale  in the ENDOR experiments for different positions of nuclear spins. At some lattice positions the off-diagonal elements of tensor $`𝐀(𝐑_k)`$, Eq. (3), are comparable with the diagonal elements, Table 1. Though the effective mass theory desribes symmetry of the hyperfine tensors incorrectly Feher , the reasonable results have been obtained in the Bloch representations Ivey1 ; Ivey2  or in Wannier representations Onffroy  of a 31P donor electron wave function.

To calculate an evolution of the electron spin density matrix with the Hamiltonian, Eq. (2), we use the approach derived by W.B. Mims Mims  for the modeling of electron spin echo envelope modulation functions (ESEEM). In a high magnetic field, we assume that the effect of the non-diagonal terms ($`S_+I_\alpha +S_{}I_\alpha `$, $`\alpha =+,,z`$) in Eq. (2) is small. These transitions are forbidden in the first order of perturbation and inversely proportional to the electron Zeeman level splitting in the second order. The truncated Hamiltonian can be written in the form

$$H_{\mathrm{trunc}}=g\mu _\mathrm{B}\mathrm{H}S_z+a_0S_zI_z^0\gamma _\mathrm{P}\mathrm{}\mathrm{H}I_z^0+\underset{k}{}\left(\underset{\alpha }{}A_{z\alpha }(𝐑_k)S_zI_\alpha ^k\gamma _{\mathrm{Si}}\mathrm{}\mathrm{H}I_z^k\right).$$
(4)
It is diagonal for the electron spin but non-diagonal for spins of 29Si nuclei. In the interaction representation with

$$H_0=g\mu _\mathrm{B}\mathrm{H}S_z+a_0S_zI_z^0\gamma _\mathrm{P}\mathrm{}\mathrm{H}I_z^0,$$
(5)
the Hamiltonian, Eq. (4), can be diagonalized for different eigenvalues $`m_s`$ of $`S_z`$ using the rotation operator

$$S=\underset{k}{}e^{i\phi _k(m_s)I_y^k}e^{i\theta _kI_z^k},$$
(6)
where $`\mathrm{tan}\theta _k=A_{zy}(𝐑_k)/A_{zx}(𝐑_k)`$, $`A_{}(𝐑_k)=\sqrt{A_{zx}^2(𝐑_k)+A_{zy}^2(𝐑_k)}`$ and $`\mathrm{tan}\phi _k(m_s)=m_sA_{}(𝐑_k)/(A_{zz}(𝐑_k)m_s\gamma _{\mathrm{Si}}\mathrm{}\mathrm{H})`$.

We consider a $`\pi /2`$ magnetic pulse as the single qubit operation, $`Op`$, and assume that initially ($`t=0^{}`$) the nuclear system is in equilibrium and unpolarized. The electron spin density matrix in the interaction picture, Eq. (5), at time t is,

$`\left(\rho _\sigma (t)\right)_{m_sm_s^{}}={\displaystyle \underset{k}{}}\{\mathrm{cos}^2\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{cos}\left({\displaystyle \frac{\omega _k(m_s^{})\omega _k(m_s)}{2}}t\right)`$
$`+\mathrm{sin}^2\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{cos}\left({\displaystyle \frac{\omega _k(m_s^{})+\omega _k(m_s)}{2}}t\right)\left\}\right(\rho _\sigma (0^+))_{m_sm_s^{}},`$ (7)
where the effective oscillation frequency of the nuclear spin located at the lattice site k with the electron spin projection $`m_s`$ is

$$\omega _k(m_s)=1/\mathrm{}\sqrt{(A_{zz}(𝐑_k)m_s\gamma _{\mathrm{Si}}\mathrm{}\mathrm{H})^2+\left(A_{}(𝐑_k)m_s\right)^2}.$$
(8)
In a high magnetic field, where

$$\mathrm{H}|A_{zz}(𝐑_k)/(2\gamma _{\mathrm{Si}}\mathrm{})|,$$
(9)
this effect can be treated as a perturbation, and Eq. (7) is simplified to

$`\left(\rho _\sigma (t)\right)_{m_sm_s^{}}={\displaystyle \underset{k}{}}\{(1{\displaystyle \frac{\left(A_{}(𝐑_k)/2\mathrm{}\right)^2}{(\gamma _{\mathrm{Si}}\mathrm{H})^2}})^2\mathrm{cos}\left({\displaystyle \frac{A_{zz}(𝐑_k)t}{2\mathrm{}}}(1{\displaystyle \frac{\left(A_{}(𝐑_k)/2\mathrm{}\right)^2}{(\gamma _{\mathrm{Si}}\mathrm{H})^2}})\right)`$
$`+{\displaystyle \frac{\left(A_{}(𝐑_k)/2\mathrm{}\right)^2}{(\gamma _{\mathrm{Si}}\mathrm{H})^2}}\mathrm{cos}\left(\gamma _{\mathrm{Si}}\mathrm{H}t(1{\displaystyle \frac{\left(A_{}(𝐑_k)/2\mathrm{}\right)^2}{(\gamma _{\mathrm{Si}}\mathrm{H})^2}})\right)\left\}\right(\rho _\sigma (0^+))_{m_sm_s^{}}.`$ (10)
The reorientation of the effective magnetic field at each nuclear spin position changes the static component of the hyperfine field at the bounded electron, the first term in Eq. (10), and admixes oscillations near the 29Si nuclear Larmor frequency, the second term in Eq. (10). The amplitude and width of the admixed band are proportional to H-2. This indicates that the effect can be reduced to any expected value by tuning of the external magnetic field. In the case of a polarized nuclear spin bath, the polarization $`p_k`$ for each nuclear spin, in analogy with a spin temperature, should be used. Therefore, Eq. (7) transforms to

$`\left(\rho _\sigma (t)\right)_{m_sm_s^{}}={\displaystyle \underset{k}{}}\left(\rho _\sigma (0^+)\right)_{m_sm_s^{}}\{\mathrm{cos}^2\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{cos}\left({\displaystyle \frac{\omega _k(m_s^{})\omega _k(m_s)}{2}}t\right)`$
$`+\mathrm{sin}^2\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{cos}\left({\displaystyle \frac{\omega _k(m_s^{})+\omega _k(m_s)}{2}}t\right)`$
$`+ip_k\mathrm{cos}\left({\displaystyle \frac{\phi _k(m_s^{})+\phi _k(m_s)}{2}}\right)\mathrm{cos}\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{sin}\left({\displaystyle \frac{\omega _k(m_s^{})\omega _k(m_s)}{2}}t\right)`$
$`ip_k\mathrm{sin}\left({\displaystyle \frac{\phi _k(m_s^{})+\phi _k(m_s)}{2}}\right)\mathrm{sin}\left({\displaystyle \frac{\phi _k(m_s^{})\phi _k(m_s)}{2}}\right)\mathrm{sin}\left({\displaystyle \frac{\omega _k(m_s^{})+\omega _k(m_s)}{2}}t\right)\}.`$ (11)
In a particular qubit implementation, the qubit density matrix will depend on the position and polarization of each 29Si spin. We estimate the value of this field reorientation effect for a 31P electron spin interacting with a single 29Si nuclear spin placed at one of the lattice sites from Table 1 in an external magnetic field, H = 1 T. The magnetic isotope 29Si, which is located in one of the four, nearest to the donor, lattice positions, experiences the largest reorientation of the effective magnetic field due to the electron spin manipulation, Fig. 1. If the nuclear spin is initially oriented parallel to $`H_{\mathrm{eff}}`$, Fig. 1(a), then after the electron spin flips, it will rotate about new axis of quantization, Fig. 1(b). To quantify the influence of this rotation on the quantum bit evolution, we calculated the norm of deviation of qubit state from ideal one Fedichkin

$$\mathrm{\Delta }\rho _\sigma (t)=\rho _\sigma (t)\rho _\sigma ^{id}(t),$$
(12)
where the deviation of the norm for single qubit is $`A=\sqrt{A_{11}A_{22}+A_{12}A_{21}}`$. As the ideal state we used the spin polarized state, Eq. (11), without high frequency components.

The simulated time evolution of the density matrix deviation, after $`Op=\pi /2`$, is shown in Fig. 2 for an initial nuclear spin polarization, $`p_k=1`$. For time scales on the order of 100 ns, the amplitude of the density matrix deviation varies from 0 to $`7\times 10^3`$, depending on the nuclear spin position.

To perform practically useful large scale quantum calculations, quantum error correction procedures should be used. The best known quantum error-correction schemes require the error rate to be less than one error per sequence of $`10^310^5`$ quantum gates Preskill ; Gott ; DiV . Our results indicate, for solid state nuclear spin quantum computers, the considered mechanism can be crucial, especially for qubits in which the 29Si isotopes are located at the nearest positions to the 31P.

For natural silicon crystal, n(29Si)$`4.76\%`$, the probability to have a 29Si in one of the four neighbors to the donor position is approximately $`18\%`$. Isotopic purification of a silicon crystal will reduce the probability to have a “wrong” qubit only. To suppress undesirable qubit dynamics, one should use high magnetic fields. It is possible to estimate the value of threshold magnetic field, using the $`\mathrm{\Delta }\rho _{\mathrm{max}}\mathrm{H}^2`$ dependence, from $`\mathrm{\Delta }\rho Q=1`$ relation, where $`Q10^4`$ is quality factor required for application of quantum error correction schemes. If we assume the effect of a whole nuclear magnetic bath is of the order of the maximal single nuclear spin effect, then the threshold value of the external magnetic field should be $`\mathrm{H}_{\mathrm{th}}9`$ T.

In the opposite case, to analyze this small effect, one should define experimental conditions in which it has considerable influence. To investigate the electron spin echo modulations, due to a nuclear spin at the E position (Table 1), the magnetic field should be, at least, of the order of H$`0.1`$ T (the electron Zeeman frequency $`\nu _Z2.8`$ GHz). Although the ESEEM, due to interaction with 29Si nuclear spins, has been observed on the irradiation-induced defects in quartz at the electron Zeeman frequency $`\nu _Z=1.7`$ GHz Kurshev , as far as we know, there was no published works on spin echo modulations for the phosphorus donor electron spin.

In the case of GaAs or InGaAs quantum dots the concentration of nuclear magnetic isotopes is $`100\%`$. The electron spin interacts with a much larger set of a nuclear spins, and we expect more complicated dynamics of system. While there are no experimental values of hyperfine constants for GaAs and InGaAs quantum dots, it may be possible to estimate them theoretically. We suppose that, in analogy with silicon, calculations in the effective mass approximation would be incorrect and one should use, at least, the LCAO model.

In conclusion, we considered the effect of the non-diagonal terms in the dipole-dipole interaction of a 31P electron spin with a system of 29Si nuclei in a high magnetic field. Electron spin manipulations generate fast oscillations of the nuclear spin bath. Due to this effect the electron spin density matrix deviates from an ideal state. In an external magnetic field of $`\mathrm{H}=1`$ T, the value of deviation is on the order of $`7\times 10^3`$. This is crucial for the implementation of practically useful quantum algorithms with electron spin qubit. The estimated threshold value of the magnetic field for a fault tolerant 31P electron spin qubit is $`\mathrm{H}_{\mathrm{th}}9`$ T.

###### Acknowledgements.

We thank S. A. Lyon, J. Nesteroff and V. Privman for helpful discussions. This research was supported by the National Security Agency and Advanced Research and Development Activity under Army Research Office contract DAAD-19-02-1-0035, and by the National Science Foundation, grants DMR-0121146 and ECS-0102500.


text:
# About rational-trigonometric deformation

## 1 Introduction

In classical mathematics there are three classes of meromorphic functions: rational, trigonometric, elliptic. According to these classes there are three types of quantum deformations: rational, trigonometric, elliptic. It turns out we can also introduce a rational-trigonometric deformation. All these deformations will be called standard.

In this short and sketched paper we discussed the simplest examples of the standard deformations in the arithmetic (number theory), geometry, differential calculus, theory functions and Lie algebras. We also introduce rational-trigonometric differential Knizhnik-Zamolodchikov (KZ) and dynamical (DD) equations. The rational-trigonometric differential KZ equations are connected with a rational-trigonometric classical $`r`$-matrix  which is a sum the simplest rational and trigonometric $`r`$-matrices depending on spectral parameter. It turns out that the rational-trigonometric differential KZ (DD) equations are sums of the rational and trigonometric differential KZ (DD) equations.

## 2 Trigonometric, rational and rational-trigonometric deformations

### 2.1 Standard deformations of numbers

Let $`z`$ be a complex number, $`z\mathrm{I}\mathrm{C}`$.

a). The trigonometric deformation (or $`q`$-deformation) of $`z`$:

$$z(z)_q=\frac{1q^z}{1q},(1)_q=1,(0)_q=0,$$
(1)
where $`q`$ is a deformation parameter, $`(z)_{q=1}=z`$. This deformation is well-known.

b). The rational deformation (or $`\eta `$-deformation) of $`z`$ :

$$z(z)_\eta =\frac{z}{1+\eta (z1)},(1)_\eta =1,(0)_\eta =0,$$
(2)
where $`\eta `$ is a deformation parameter, $`(z)_{\eta =0}=z`$.

c). The rational-trigonometric deformation (or $`(q,\eta )`$-deformation) of $`z`$ :

$$z(z)_{q\eta }=\frac{(z)_q}{1+\eta (z1)_q},(1)_{q\eta }=1,(0)_{q\eta }=0,$$
(3)
where $`q`$ and $`\eta `$ are deformation parameters, $`(z)_{q=1,\eta }=(z)_\eta `$, $`(z)_{q,\eta =0}=(z)_q`$.

Remark. There is an elliptic deformation of $`z`$ (for example, see )

### 2.2 Standard deformations of two-dimensional plane

Let $`x`$ and $`y`$ be two commuting coordinate variables, i.e.

$$xyyx=0(\mathrm{the}(x,y)\mathrm{plane}).$$
(4)
a). The trigonometric deformation of the $`(x,y)`$-plane:

$$xyqyx=0(\mathrm{the}\mathrm{Manin}^{}\mathrm{s}\mathrm{plane}).$$
(5)
b). The rational deformation of the $`(x,y)`$-plane:

$$xyyx=\eta y^2.$$
(6)
c). The rational-trigonometric deformation of the $`(x,y)`$-plane:

$$xyqyx=\eta y^2.$$
(7)
All these deformation of the $`(x,y)`$-plane are well-known.

Remark. No elliptic deformation of the $`(x,y)`$-plane is known. It is an open problem.

### 2.3 Standard deformations of the differential calculus

Let $`_x:=\frac{}{x}`$ be the usual derivative.

a). The trigonometric deformation of $`_x`$:

$$_x^{(q)}f(x)_q=\frac{f(qx)f(x)}{qxx}=\frac{f(qx)f(x)}{x(q1)}.$$
(8)
b). The rational deformation of $`_x`$:

$$_x^{(\eta )}f(x)_q=\frac{f(x+\eta )f(x)}{x+\eta x}=\frac{f(x+\eta )f(x)}{\eta }.$$
(9)
The deformations a) and b) for $`_x`$ are well-known.

c). The rational-trigonometric deformation of $`_x`$:

$$_x^{(q\eta )}f(x)_q=\frac{f(qx+\eta )f(x)}{qx+\eta x}=\frac{f(qx+\eta )f(x)}{x(q1)+\eta }.$$
(10)
All these deformations of $`_x`$ have the following properties:

(i) they are defined by the general formulas

$$_x^{(df)}f(x)=\frac{\mathrm{}^{(df)}f(x)}{\mathrm{}^{(df)}x}=\frac{f(x^{})f(x)}{x^{}x},$$
(11)
(ii) they satisfy the Leibniz rule

$$_x^{(df)}f(x)\phi (x)=\left(_x^{(df)}f(x)\right)\phi (x)+f(x^{})\left(_x^{(df)}\phi (x)\right),$$
(12)
where $`x^{}=qx`$ for the $`q`$-deformation, $`x^{}=x+\eta `$ for the $`\eta `$-deformation and $`x^{}=qx+\eta `$ for the $`(q,\eta )`$-deformation.

Remark. There is an explicite formula for the elliptic deformation of $`_x`$.

### 2.4 Standard deformations of functions

Let $`\mathrm{exp}(x)`$ be the usual exponential of $`x`$:

$$\mathrm{exp}(x)=1+x+\frac{1}{2!}x^2+\mathrm{}+\frac{1}{n!}x^n+\mathrm{}.$$
(13)
and let $`F_{n,m}(x)`$ be a standard hypergeometric series:

$$F_{n,m}(x)=\underset{k=0}{\overset{\mathrm{}}{}}\frac{(a_1)_{(k)}(a_2)_{(k)}\mathrm{}(a_n)_{(k)}}{k!(b_1)_{(k)}(b_2)_{(k)}\mathrm{}(b_m)_{(k)}}x^k,$$
(14)
where

$$(a)_{(k)}=a(a+1)\mathrm{}(a+k1).$$
(15)
a). The trigonometric deformation of $`\mathrm{exp}(x)`$ and $`F_{n,m}(x)`$.

$$\mathrm{exp}_q(x)=1+x+\frac{1}{(2)_q!}x^2+\mathrm{}+\frac{1}{(n)_q!}x^n+\mathrm{},$$
(16)
where

$$(n)_q!=(1)_q(2)_q\mathrm{}(n)_q.$$
(17)
If we replace the parameters $`a_i`$, $`b_j`$, and $`k!`$ in (14) by the $`q`$-analogs $`(a_i)_q`$, $`(b_j)_q`$ and $`(k)_q!`$ we obtain the basic hypergeometric ($`q`$-hypergeometric) series $`F_{n,m}^{(q)}(z)`$.

b). The rational deformation of $`\mathrm{exp}(x)`$ and $`F_{n,m}(x)`$.

$$\mathrm{exp}_\eta (x)=1+x+\frac{1}{(2)_\eta !}x^2+\mathrm{}+\frac{1}{(n)_\eta !}x^n+\mathrm{},$$
(18)
where

$$(n)_\eta !=(1)_\eta (2)_\eta \mathrm{}(n)_\eta .$$
(19)
If we replace the parameters $`a_i`$, $`b_j`$, and $`k!`$ in (14) by the $`\eta `$-analogs $`(a_i)_\eta `$, $`(b_j)_\eta `$ and $`(k)_\eta !`$ we obtain the basic hypergeometric ($`\eta `$-hypergeometric) series $`F_{n,m}^{(\eta )}(z)`$.

c). The rational-trigonometric deformation of $`\mathrm{exp}(x)`$ and $`F_{n,m}(x)`$.

$$\mathrm{exp}_{q\eta }(x)=1+x+\frac{1}{(2)_{q\eta }!}x^2+\mathrm{}+\frac{1}{(n)_{q\eta }!}x^n+\mathrm{},$$
(20)
where

$$(n)_{q\eta }!=(1)_{q\eta }(2)_{q\eta }\mathrm{}(n)_{q\eta }.$$
(21)
The replacement of $`a_i`$, $`b_j`$, and $`k!`$ in (14) by the $`(q,\eta )`$-analogs $`(a_i)_{q\eta }`$, $`(b_j)_{q\eta }`$ and $`(k)_{q\eta }!`$ gives us the $`(q,\eta )`$-hypergeometric series

$$F_{n,m}^{(q,\eta )}(x)=\underset{k=0}{\overset{\mathrm{}}{}}\frac{((a_1)_{q\eta })_{(k)}((a_2)_{q\eta })_{(k)}\mathrm{}((a_n)_{q\eta })_{(k)}}{(k)_{q\eta }!((b_1)_{q\eta })_{(k)}((b_2)_{q\eta })_{(k)}\mathrm{}((b_m)_{q\eta })_{(k)}}x^k,$$
(22)
where

$$((a)_{q\eta })_{(k)}=(a)_{q\eta }(a+1)_{q\eta }\mathrm{}(a+k1)_{q\eta }.$$
(23)
Setting here $`q=1`$ we obtain the $`\eta `$-hypergeometric series $`F_{n,m}^{(\eta )}(x)`$. We can also introduce the $`\eta `$\- and ($`q,\eta `$)- analogs of other special functions.

All deformed exponentials $`\mathrm{exp}_q(x)`$, $`\mathrm{exp}_\eta (x)`$ and $`\mathrm{exp}_{q\eta }(x)`$ can be obtained from the functional equation $`f_1(x+y)=f_2(y)f_3(x)`$ with regular functions $`f_i(z)`$ satisfying the initial conditions $`f_i(0)=1`$, where the variables $`x`$ and $`y`$ satisfy the relations (5)–(7) (see ).

Remark. Analogous formulas exist also for the elliptic case they can be found in .

### 2.5 Standard deformations of universal enveloping algebras

Let $`g`$ be a finite-dimensional Lie algebra, $`g[u,u^1]`$ be a loop algebra and $`g[u]`$ be a non-negative loop algebra over $`g`$. We denote by $`U(g)`$, $`U(g[u])`$ and $`U(g[u,u^1])`$ their universal enveloping algebras.

a). The trigonometric deformation ($`q`$-deformation) of $`U(g)`$, $`U(g[u])`$ and $`U(g[u,u^1])`$ are well-known. They are denoted by $`U_q(g)`$, $`U_q(g[u])`$ and $`U_q(g[u,u^1])`$ .

b). The rational deformation of $`U_q(g[u])`$ and $`U_q(g[u,u^1])`$ are well-known. They are Yangian $`Y_\eta (g)`$  and its double $`DY_\eta (g)`$ .

c). The rational-trigonometric deformation of $`U_q(g[u])`$ and $`U_q(g[u,u^1])`$ are Drinfeldian $`D_{q\eta }(g)`$ ,  and its double $`DD_{q\eta }(g)`$.

Remark. An elliptic deformation exists only for $`U(sl_n[u])`$. It is well-known too.

## 3 Rational-trigonometric differential Knizhnik-Zamolodchikov and dynamical equations

We consider the case $`g=gl_M`$ although results of this section are also valid for an arbitrary simple complex Lie algebra $`g`$. Let $`e_{ab}`$, $`a,b=1,\mathrm{},M`$, be the standard Cartan-Weyl basis of $`gl_M`$: $`[e_{ab},e_{cd}]=\delta _{bc}e_{ad}\delta _{ad}e_{cb}`$. The element $`C_2:=\underset{a,b=1}{\overset{N}{}}e_{ab}e_{ba}U(gl_M)`$ is a $`gl_M`$-scalar, i.e. $`[C_2,x]=0`$ for any $`xgl_M`$, and it is called the second order Casimir element. The element $`\mathrm{\Omega }:=\frac{1}{2}\left(\mathrm{\Delta }(C_2)C_2\mathrm{id}\mathrm{id}C_2\right)=\underset{a,b=1}{\overset{M}{}}e_{ab}e_{ba}U(gl_M)U(gl_M)`$, where $`\mathrm{\Delta }`$ is a trivial co-product $`\mathrm{\Delta }(x)=x\mathrm{id}+\mathrm{id}x`$ ($`xgl_M`$), is called the Casimir two-tensor. The two-tensor can be presented in the form $`\mathrm{\Omega }=`$ $`\mathrm{\Omega }^++\mathrm{\Omega }^{}`$, where $`\mathrm{\Omega }^+=\frac{1}{2}\underset{a}{}e_{aa}e_{aa}+\underset{1a<bM}{}e_{ab}e_{ba}`$ and $`\mathrm{\Omega }^{}=\frac{1}{2}\underset{a}{}e_{aa}e_{aa}+\underset{1a<bM}{}e_{ba}e_{ab}`$. Note that $`(\omega \omega )(\mathrm{\Omega }^\pm )=\mathrm{\Omega }^{}`$, where $`\omega `$ is the Cartan automorphism: $`\omega (e_{ab})=e_{ba}`$.

For any $`xU(gl_M)`$ we set $`x_{(i)}=\underset{Ntimes}{\overset{ith}{\underset{}{\mathrm{id}\mathrm{}\text{id}𝑥\text{id}\mathrm{}\text{id}}}}`$. We consider $`U(gl_M)`$ as a subalgebra of $`(U(gl_M))^N`$, the embedding $`U(gl_M)(U(gl_M))^N`$ being given by the $`N`$-fold co-product, that is $`x\mathrm{\Delta }^N(x)=x_{(1)}+\mathrm{}+x_{(N)}`$ for any $`xgl_M`$.

For a nonzero complex number $`\kappa `$ we consider differential operators $`_{z_1}^{(r)},\mathrm{},_{z_N}^{(r)}`$ and $`_{z_1}^{(t)},\mathrm{},_{z_N}^{(t)}`$ with coefficient in $`(U(gl_M))^N`$ depending on complex variable $`z_1,\mathrm{},z_N`$ and $`\lambda _1,\mathrm{},\lambda _M`$ (see  and ):

$$_{z_i}^{(r)}(z;\lambda )=\kappa \frac{}{z_i}\underset{a=1}{\overset{M}{}}\lambda _a(e_{aa})_{(i)}\underset{\genfrac{}{}{0pt}{}{j=1}{ji}}{\overset{N}{}}\frac{\mathrm{\Omega }_{(ij)}}{z_iz_j},$$
(24)
$`_{z_i}^{(t)}(z;\lambda )`$ $`=`$ $`\kappa z_i{\displaystyle \frac{}{z_i}}{\displaystyle \underset{a=1}{\overset{M}{}}}(\lambda _ae_{aa})(e_{aa})_{(i)}{\displaystyle \underset{\genfrac{}{}{0pt}{}{j=1}{ii}}{\overset{N}{}}}{\displaystyle \frac{z_i\mathrm{\Omega }_{(ij)}^++z_j\mathrm{\Omega }_{(ij)}^{}}{z_iz_j}}.`$ (25)
The operators $`_{z_1}^{(r)},\mathrm{},_{z_N}^{(r)}`$ (resp.$`_{z_1}^{(t)},\mathrm{},_{z_N}^{(t)}`$) are called the rational (resp. trigonometric) Knizhnik-Zamolodchikov (KZ) operators. We set

$$_{z_i}^{(rt)}(z;\lambda )=\mathrm{}_{z_i+\frac{\eta }{\mathrm{}}}^{(t)}(z+\frac{\eta }{\mathrm{}};(1+\frac{\eta }{\mathrm{}})\lambda ).$$
(26)
It is easy to see that

$$_{z_i}^{(rt)}(z;\lambda )=\mathrm{}_{z_i}^{(t)}(z;\lambda )+\eta _{z_i}^{(r)}(z;\lambda ).$$
(27)
The rational (trigonometric or rational-trigonometric) KZ equations is the system of the differential equations

$$_{z_i}^{(\alpha )}(z;\lambda )u(z;\lambda )=0,i=1,\mathrm{},N;\alpha =(r)((t)\mathrm{or}(rt))$$
(28)
for a function $`u(z;\lambda ):=u(z_1,\mathrm{},z_N;\lambda _1,\mathrm{},\lambda _M)`$ taking values in an $`N`$-fold tensor product of $`gl_M`$-modules.

We also consider the differential operators $`D_{\lambda _1}^{(r)},\mathrm{},D_{\lambda _M}^{(r)}`$ and $`D_{\lambda _1}^{(t)},\mathrm{},D_{\lambda _M}^{(t)}`$ with coefficient in $`(U(gl_M))^N`$ depending on complex variables $`z_1,\mathrm{},z_N`$ and $`\lambda _1,\mathrm{},\lambda _M`$ :

$$D_{\lambda _a}^{(r)}(z;\lambda )=\kappa \frac{}{\lambda _a}\underset{i=1}{\overset{N}{}}z_i(e_{aa})_{(i)}\underset{\genfrac{}{}{0pt}{}{b=1}{ba}}{\overset{M}{}}\frac{e_{ab}e_{ba}e_{aa}}{\lambda _a\lambda _b},$$
(29)
$`D_{\lambda _a}^{(t)}(z;\lambda )`$ $`=`$ $`\kappa \lambda _a{\displaystyle \frac{}{\lambda _a}}+{\displaystyle \frac{a_{aa}^2}{2}}{\displaystyle \underset{i=1}{\overset{N}{}}}z_i(e_{aa})_i`$ (30)
$`{\displaystyle \underset{b=1}{\overset{M}{}}}{\displaystyle \underset{1i<jN}{}}(e_{ab})_{(i)}(e_{ba})_{(j)}{\displaystyle \underset{\genfrac{}{}{0pt}{}{b=1}{ba}}{\overset{M}{}}}{\displaystyle \frac{\lambda _b(e_{ab}e_{ba}e_{aa})}{\lambda _a\lambda _b}}.`$
Remind that $`e_{ab}=\underset{i=1}{\overset{M}{}}(e_{ab})_{(i)}`$. The operators $`D_{z_1}^{(r)},\mathrm{}D_{z_M}^{(r)}`$ (resp. $`D_{z_1}^{(t)},\mathrm{}D_{z_M}^{(t)}`$) are called the rational (resp. trigonometric) differential dynamic (DD) operators. We set

$$D_{\lambda _a}^{(rt)}(z;\lambda ):=\mathrm{}D_{\lambda _a+\frac{\eta }{\mathrm{}}}^{(t)}((1+\frac{\eta }{\mathrm{}})z;\lambda +\frac{\eta }{\mathrm{}}).$$
(31)
It is easy to see that

$$D_{z_i}^{(rt)}(z;\lambda )=\mathrm{}D_{\lambda _a}^{(t)}(z;\lambda )+\eta D_{\lambda _a}^{(r)}(z;\lambda ).$$
(32)
The rational (trigonometric or rational-trigonometric) DD equations is the system of the differential equations

$$D_{\lambda _a}^{(\alpha )}(z;\lambda )u(z;\lambda )=0,a=1,\mathrm{},M;\alpha =(r)((t)\mathrm{or}(rt))$$
(33)
for a function $`u(z;\lambda ):=u(z_1,\mathrm{},z_N;\lambda _1,\mathrm{},\lambda _M)`$ taking values in an $`N`$-fold tensor product of $`gl_M`$-modules.

¿From (27), (32) and Theorem 5.8 of the paper  for any $`i=1,\mathrm{},N`$ and $`a=1,\mathrm{},M`$ we have the duality:

$$_{z_i}^{(rt)}(z;\lambda )_MD_{z_i}^{(rt)}(\lambda ;z)_N,D_{\lambda _a}^{(rt)}(z;\lambda )_M_{\lambda _a}^{(rt)}(\lambda ;z)_N,$$
(34)
where the index $`M`$ ($`N`$) is connected with the Lie algebra $`gl_M`$($`gl_N`$).

Acknowledgment. This work was supported by Russian Foundation for Fundamental Research, grant No. RFBR-02-01-00668, and CRDF RMI-2334-MO-02.


text:
# 1 Introduction

## 1 Introduction

The massive Schwinger model is QED in 1+1 dimensions with one massive fermion species . It is, one the one hand, the simplest non-trivial gauge field theory. On the other hand, it nevertheless shares some properties with more complicated gauge field theories like, e.g., QCD in 3+1 dimensions. The model has, for instance, a non-trivial vacuum , and the fermions are confined in the sense that they do not show up as physical particles in the particle spectrum .

The Lagrangian density of the model is

$$L=\overline{\mathrm{\Psi }}(i/gA/m)\mathrm{\Psi }\frac{1}{4}F_{\mu \nu }F^{\mu \nu }.$$
(1)
Here, both the fermion mass $`m`$ and the coupling constant $`g`$ have the dimension of a mass, therefore the dimensionless parameter $`m/g`$ effectively controls the behaviour of the theory.

One interesting feature of the model is the formation of bound states. For small coupling (i.e., large $`m/g`$), the bound states are just fermion–antifermion bound states that are weakly bound by the 1+1 dimensional Couloumb potential . Therefore, the masses of the bound states are $`2m`$ plus some small (positive) binding energies and should approach $`2m`$ in the limit of vanishing coupling $`g`$.

For large coupling, the model is best described in an equivalent bosonic formulation . In fact, for a massless fermion ($`m=0`$), the model may be solved exactly, and it is equivalent to the theory of one free, massive boson field (“Schwinger boson”) with Schwinger mass $`\mu _0=g/\sqrt{\pi }`$ . For small but non-zero fermion mass ($`m0`$), this Schwinger boson turns into an interacting particle, and its mass acquires corrections which can be calculated perturbatively in $`m`$ . This Schwinger boson may be identified with the lowest bound state, the vector state. Further, a bound state of two Schwinger bosons is formed within the bosonic formulation of the theory , and this bound state of two bosons may be identified with the next lowest bound state, the scalar state.

These two lowest bound states remain stable for all values of the parameter $`m/g`$. In general, all higher bound states become unstable for sufficiently large coupling (small $`m/g`$) . For vanishing vacuum angle ($`\theta =0`$) there exists a third stable bound state , which may be interpreted as a bound state of one Schwinger boson and one two-boson bound state within a bosonic formulation . Here we will consider the case of vanishing vacuum angle (because all lattice data are for this case), but we will be concerned only with the two lowest bound states, because to our knowledge there are no lattice data available for the third bound state.

For the calculation we will use perturbative expansions of the bosonic formulation of the theory for different choices of the normal-ordering mass, together with some consistency conditions among these different perturbation expansions. Qualitatively this has already been discussed in , therefore we will describe the calculation only briefly before presenting explicit results.

## 2 Vector mass calculation

The theory of Eq. (1) is equivalent to the bosonic theory  (for conventions we refer to )

$$L=N_{\mu _0}\left[\frac{1}{2}\mathrm{\Phi }(\mathrm{}\mu _0^2)\mathrm{\Phi }+\frac{e^\gamma }{2\pi }\mu _0m\left(\mathrm{cos}(\sqrt{4\pi }\mathrm{\Phi })1\right)\right].$$
(2)
Here $`\mu _0=g/\pi ^{1/2}`$ is the Schwinger mass of the massless ($`m=0`$) Schwinger model, $`\gamma =0.5772`$ is the Euler constant and $`N_\mu `$ denotes normal-ordering w.r.t. $`\mu `$. Further, a constant vacuum energy density has been subtracted, as usual.

By making use of the well-known normal-ordering relation

$$N_{\overline{\mu }}e^{\pm i\beta \mathrm{\Phi }(x)}=(\frac{\mu }{\overline{\mu }})^{\frac{\beta ^2}{4\pi }}N_\mu e^{\pm i\beta \mathrm{\Phi }(x)}$$
(3)
(see ), we may rewrite the Lagrangian density (2) for arbitrary normal-ordering mass $`\overline{\mu }`$ like

$$L=N_{\overline{\mu }}\left[\frac{1}{2}\mathrm{\Phi }(\mathrm{}\mu ^2)\mathrm{\Phi }+\frac{e^\gamma }{2\pi }\overline{\mu }m\left(\mathrm{cos}(\sqrt{4\pi }\mathrm{\Phi })1\right)\right].$$
(4)
By shifting the quadratic (in $`\mathrm{\Phi }`$) part from the interaction part $`\mathrm{cos}\sqrt{4\pi }\mathrm{\Phi }`$ to the free part of the Lagrangian (4), and by insisting that the normal-ordering mass is equal to the total mass of the free Lagrangian, we arrive at

$$L=N_{\overline{\mu }_0}\left[\frac{1}{2}\mathrm{\Phi }(\mathrm{}\overline{\mu }_0^2)\mathrm{\Phi }+\frac{e^\gamma }{2\pi }m\overline{\mu }_0\left(\mathrm{cos}(\sqrt{4\pi }\mathrm{\Phi })+2\pi \mathrm{\Phi }^21\right)\right]$$
(5)
with

$$\overline{\mu }_0=e^\gamma m+\sqrt{\mu _0^2+e^{2\gamma }m^2}.$$
(6)
The resulting Lagrangian (5) has a slightly more complicated interaction term (containing a $`\mathrm{\Phi }^2`$ term in addition to the $`\mathrm{cos}\sqrt{4\pi }\mathrm{\Phi }`$ term), and the dimensionless expansion parameter is now

$$\lambda \frac{m}{\overline{\mu }_0}.$$
(7)
This new expansion parameter $`\lambda `$ has the advantage that it tends to a finite and rather small value even for large $`m/g`$.

The lowest bound-state mass $`M_1`$ (the “vector” state mass or Schwinger mass) may be calculated from the two-point function $`\mathrm{\Phi }(x)\mathrm{\Phi }(y)`$ for both Lagrangians (2) and (5) by a perturbation expansion in the respective interaction terms. Via a Fourier transformation the mass pole may be determined. For ordinary mass perturbation theory (i.e., from the Lagrangian (2)), the result up to second order is

$$(M_1^{m,2})^2=\mu _0^2+2e^\gamma \mu _0m+e^{2\gamma }m^2(A_++A_{}),$$
(8)
where

$$A_\pm =_0^{\mathrm{}}𝑑rr\left[e^{2K_0(r)}1+I_0(r)(\pm e^{2K_0(r)}1+2K_0(r))\right]$$
(9)
$$A_+=0.6599,A_{}=1.7277$$
(10)
see e.g. . The notation $`M_1^{m,2}`$ means that the subscript $`1`$ stands for the first (lowest) bound state mass, and the superscript $`m,2`$ stands for the $`m`$ perturbation expansion up to second order. Here, $`\mu _0^2`$ is the order zero result, and there are both first and second order contributions. This result describes the vector state mass quite well for small fermion mass $`m`$, see Fig. 1, but it diverges, of course, in the limit of very large $`m`$.

Starting from the Lagrangian (5), the order zero result is $`\overline{\mu }_0`$. There is no first order contribution (because there is no $`\mathrm{\Phi }^2`$ term in the interaction Lagrangian), and the total mass up to second order is

$$(M_1^{\lambda ,2})^2=\overline{\mu }_0^2+e^{2\gamma }m^2(A_++A_{})2e^{2\gamma }m^2$$

$$=\mu ^2+2e^\gamma \mathrm{cos}\theta m\sqrt{\mu ^2+e^{2\gamma }\mathrm{cos}^2\theta m^2}+e^{2\gamma }m^2(A_+\mathrm{cos}2\theta +A_{}),$$
(11)
see . Please observe that one piece of the second order (in $`\lambda `$) contribution, namely the term $`2e^{2\gamma }m^2`$, has precisely cancelled the $`o(m^2)`$ piece of the lowest order (in $`\lambda `$) contribution, $`\overline{\mu }_0^2`$. This cancellation ensures that $`(M_1^{\lambda ,2})^2`$, when Taylor expanded in $`m`$ up to second order, coincides with the second order result of mass perturbation theory, $`(M_1^{m,2})^2`$. $`(M_1^{\lambda ,2}/m)`$ reaches a finite value (which is about 55% above the true value $`2`$) in the limit of large $`m`$. Further, $`M_1^{\lambda ,2}`$ describes the data quite well for small $`m/g`$, but already for $`m/g=1`$ it is about 30% above the true value. Therefore, a further improvement would be desireable. This aim can be achieved by requiring that the result of the $`\lambda `$ perturbation expansion, when Taylor-expanded in $`m`$, should agree with the result from the $`m`$ perturbation expansion also in higher orders. This requirement is perfectly sensible, because $`\lambda `$ has a well-defined Taylor expansion into $`m`$ to all orders, and the physical results from both formulations of the theory should agree after all. A further Taylor expansion in $`m`$ of the mass $`(M_1^{\lambda ,2})^2`$ contains, e.g., the term

$$\overline{\mu }_0^2|_{m^3}=e^{3\gamma }\frac{m^3}{\mu _0},$$
(12)
in third order in $`m`$, and infinitely many higher order terms. Now the crucial point is that such terms cannot be produced by the perturbation theory in $`m`$, because Wick contractions of the interaction term $`\mathrm{cos}\sqrt{4\pi }\mathrm{\Phi }`$ all must contain integrals of a special type of transcendental functions (the $`E_\pm `$ functions of ). Therefore, the higher orders in the $`\lambda `$ perturbation expansion must contain terms that precisely cancel terms like (12) when Taylor-expanded in $`m`$, as we have already observed for the second order contribution. The term (12) is, e.g., precisely cancelled by the third order term (in $`\lambda `$ perturbation expansion)

$$\delta (M_1^{\lambda ,3})^2=e^{3\gamma }\overline{\mu }_0^2\lambda ^3=e^{3\gamma }\frac{m^3}{\overline{\mu }}$$
(13)
which must, therefore, be present in the $`\lambda `$ perturbation expansion. This cancellation has to continue in higher orders, and the Taylor coefficients of higher order terms, like (13), have to be cancelled, too. This cancellation condition enables us to perform a further partial higher order calculation in $`\lambda `$ perturbation expansion up to some order $`N`$ by employing the equation

$$(M_1^{\lambda ,p})^2(M_1^{\lambda ,2})^2\overline{\mu }_0^2\underset{l=3}{\overset{N}{}}\frac{e^{l\gamma }}{l!}c_l\lambda ^l\stackrel{!}{=}(M_1^{m,2})^2+O(m^{N+1})$$
(14)
where a Taylor expansion in $`m`$ up to order $`N`$ is understood on the l.h.s. The coefficients $`c_l`$ can be determined recursively, which can be easily done with the help of Mathematica. We show the result for $`N=30`$ in Fig. 1 and in Table 1. As can be seen, the improvement for intermediate values of $`m/g`$ is dramatic. For $`m/g=1`$, e.g., the result deviates from the lattice result by less than 0.5%. In the limit of very large $`m`$, the deviation from the exact result $`2m`$ is about 8%.

Within $`\lambda `$ perturbation theory, the higher order terms that we have found from condition (14) correspond to disconnected products of Wick contractions, where in each connected component one $`\mathrm{cos}\mathrm{\Phi }`$ term is contracted with one or several $`\mathrm{\Phi }^2`$ terms of the interaction Lagrangian. There exist other types of terms which do not appear in the $`m`$ perturbation expansion (contractions of two or more $`\mathrm{cos}\mathrm{\Phi }`$ with at least one $`\mathrm{\Phi }^2`$), but these may be detected in a way similar to the above argument only when compared with higher orders in $`m`$ perturbation theory. Unfortunately, already a full third order calculation (both in $`m`$ and in $`\lambda `$ perturbation theory) is rather involved.

## 3 Scalar mass calculation

The scalar mass may be calculated from the four-point function $`_{l=1}^4\mathrm{\Phi }(x_l)`$ by focusing on the $`s`$-channel contribution and by again determining the mass pole. Here, we want to use the results of the last section for the vector mass, i.e., we want to use the masses $`(M_1^{m,2})^2`$ and $`(M_1^{\lambda ,p})^2`$ for the respective perturbation theories, because they are better approximations to the physical masses than the bare masses $`\mu _0`$ resp. $`\overline{\mu }_0`$. However, naively using these masses within the Lagrangians (2) and (5) would lead to the problem that the propagating masses are no longer equal to the normal-ordering masses, which would produce tadpole contributions. Therefore, we change the normal-ordering masses in both Lagrangians and use for the subsequent calculations the Lagrangians

$$L=N_{M_1}\left[\frac{1}{2}\mathrm{\Phi }(\mathrm{}M_1^2)\mathrm{\Phi }+\frac{e^\gamma }{2\pi }M_1m\left(\mathrm{cos}(\sqrt{4\pi }\mathrm{\Phi })1\right)\right].$$
(15)
and

$$L=N_{M_1}\left[\frac{1}{2}\mathrm{\Phi }(\mathrm{}M_1^2)\mathrm{\Phi }+\frac{e^\gamma }{2\pi }mM_1\left(\mathrm{cos}(\sqrt{4\pi }\mathrm{\Phi })+2\pi \mathrm{\Phi }^21\right)+e^\gamma m(\overline{\mu }_0M_1)\mathrm{\Phi }^2\right]$$
(16)
where $`M_1`$ will be finally chosen equal to $`M_1^{m,2}`$ in the first case and equal to $`M_1^{\lambda ,p}`$ in the second case. Please observe that the $`\mathrm{\Phi }^2`$ part of the interaction Lagrangian in (16) is not changed by the renormal-ordering and that we do not introduce a further shift such that all $`\mathrm{\Phi }^2`$ contributions belong to the free part of the Lagrangian. Here we just remark that this choice is perfectly legitimate, i.e., the resulting perturbation theory is as well-defined as the others. Later on it will become clear why this choice is even advisable for our purposes.

Now, we have to calculate the bosonic four-point function, perform the Fourier transformation, focus on the $`s`$-channel contribution and find the mass pole in the kinematic variable $`s`$. The calculation is essentially equivalent to the calculations within mass perturbation theory that have been done in , therefore we just quote the result. One finally arrives at the pole mass equation (or gap equation) for the scalar mass $`M_2`$ in leading order

$$1=\frac{e^\gamma m}{M_1}\frac{4}{y^2}\frac{1}{\sqrt{\frac{4}{y^2}1}}\mathrm{arctan}\frac{1}{\sqrt{\frac{4}{y^2}1}}$$
(17)
where $`y=M_2/M_1`$. It turns out that in leading order the gap equation is formally the same for both perturbation theories, the only difference being the different choices for $`M_1`$. The $`M_1`$ masses have been determined fully up to second order in the respective perturbation expansions. It turns out that the above gap equation induces a leading order contribution to the mass ratio $`y`$ which is of second order (see Eq. (19) below), therefore it is consistent to use the gap equation in leading order only. Further, we understand now why we had to choose the Lagrangian (16) as it stands. When $`M_1=M_1^{m,2}`$ in (further renormal-ordered) $`m`$ perturbation theory and $`M_1=M_1^{\lambda ,p}`$ in (further renormal-ordered) $`\lambda `$ perturbation theory, then it is automatically guaranteed that the resulting masses $`M_2`$ will agree up to order $`N`$ when Taylor-expanded in $`m`$. Here we should be somewhat clearer about the Taylor expansion. The point is that we want to use the gap equation (17) directly for a numerical evaluation, which is then no longer accessible to a Taylor expansion. However, it is easy to generate an expansion in the coupling constant from the gap equation (17). Just rewrite it as

$$ϵ=4\frac{e^\gamma m}{M_1}\frac{1}{\sqrt{4ϵ^2}}\mathrm{arccot}\frac{ϵ}{\sqrt{4ϵ^2}},ϵ^24y^2$$
(18)
and solve this equation iteratively, starting with $`ϵ=0`$ on the right-hand side. For the resulting expressions, we demand that they must agree up to order $`N`$ when Taylor-expanded in $`m`$.

There are several reasons for using the gap equation (17) directly for a numerical evaluation instead of a further expansion (iteration of Eq. (18)). Firstly, the leading second order contribution in $`m`$ has a rather large coefficient ,

$$\frac{M_2}{M_1}2+24.625\frac{m^2}{g^2}$$
(19)
and, therefore, approximates the lattice results rather badly already for quite small values of $`m`$. Secondly, after a further expansion, it is not even guaranteed that $`y<2`$ (as it has to be for a bound state), whereas $`y<2`$ holds automatically for the gap equation (17).

Remark: Please observe that for an iteration of Eq. (18) in leading order, as has been done in , the further normal-ordering of Eq. (15) is unnecessary, because it only changes results in higher orders.

The gap equation has been evaluated numerically with the help of Mathematica. Already within mass perturbation theory (i.e., for $`M_1=M_1^{m,2}`$) the results of the numerical solution are quite reasonable up to $`m/g1/2`$ and are much better than the leading order result (19). For the Lagrangian (16), i.e., for $`M_1=M_1^{\lambda ,p}`$, the results are shown in Table 1, again for the choice $`N=30`$. The results agree very well with the lattice data for small $`m`$. For $`m/g1`$, the discrepancy is about 4%, and in the limit of very large $`m`$ the deviation from the exact result $`M_2=2m`$ is about 9%.

## 4 Summary

In this paper we have performed an improved calculation of the two lowest bound-state masses of the massive Schwinger model (the vector and the scalar bound-state mass). For this purpose we used the bosonized version of the theory with two different normal-ordering prescriptions and with two slightly different interaction terms in the Lagrangian, see Eqs. (2) and (5). In addition, we used a consistency condition between the two perturbation theories in order to partially determine higher order contributins to the second ($`\lambda `$) perturbation theory, see Eq. (14). The resulting masses agree very well with the results from recent lattice calculations  for small and intermediate values of the parameter $`m/g`$. In the limit of very large $`m/g`$ they deviate from the known exact results by less than 10%. There are additional results on the masses $`M_1`$ and $`M_2`$ in the literature, but we chose to compare with the lattice results of  as they seem to be the most accurate numerical estimates which are available up to now. For additional results on these masses with various methods, please consult  and the literature cited there (some further references can also be found in ).

In principle, the results obtained in this paper could be further improved by methods that are analogous to those used in this paper. However, any further improvement requires a full third order calculation of the vector mass $`M_1`$ within mass perturbation theory, which is quite difficult to do, as already mentioned. Further, the methods emlployed in this paper may, in principle, be used for the calculation of further observables in the massive Schwinger model, like, e.g., higher masses or decay widths and scattering cross sections. Also, it should be possible to use a variant of this renormal-ordered perturbation theory for a discussion of the multi-flavour Schwinger model where ordinary mass perturbation theory fails due to infrared divergencies .

Acknowledgment: The author acknowledges support from the Austrian START award project FWF-Y-137-TEC of N.J. Mauser. Further thanks are due to M. Frank for help with the Mathematica programming.