format: text

text:
# Searching for Supersymmetric Dark Matter11footnote 1Summary of invited talk at CosPA 2002, 2002 International Symposium on Cosmology and Particle Astrophysics, May31 - June 2, 2002, National Taiwan University Taipei, Taiwan.

## 1 Introduction

Although there are many reasons for considering supersymmetry as a candidate extension to the standard model of strong, weak and electromagnetic interactions, one of the most compelling is its role in understanding the hierarchy problem namely, why/how is $`m_WM_P`$. One might think naively that it would be sufficient to set $`m_WM_P`$ by hand. However, radiative corrections tend to destroy this hierarchy. For example, one-loop diagrams generate

$$\delta m_W^2=𝒪\left(\frac{\alpha }{\pi }\right)\mathrm{\Lambda }^2m_W^2$$
(1)
where $`\mathrm{\Lambda }`$ is a cut-off representing the appearance of new physics, and the inequality in (1) applies if $`\mathrm{\Lambda }10^3`$ TeV, and even more so if $`\mathrm{\Lambda }m_{GUT}10^{16}`$ GeV or $`M_P10^{19}`$ GeV. If the radiative corrections to a physical quantity are much larger than its measured values, obtaining the latter requires strong cancellations, which in general require fine tuning of the bare input parameters. However, the necessary cancellations are natural in supersymmetry, where one has equal numbers of bosons $`B`$ and fermions $`F`$ with equal couplings, so that (1) is replaced by

$$\delta m_W^2=𝒪\left(\frac{\alpha }{\pi }\right)|m_B^2m_F^2|.$$
(2)
The residual radiative correction is naturally small if $`|m_B^2m_F^2|<1\mathrm{TeV}^2`$.

In order to justify the absence of interactions which can be responsible for extremely rapid proton decay, it is common in the minimal supersymmetric standard model (MSSM) to assume the conservation of R-parity. If R-parity, which distinguishes between “normal” matter and the supersymmetric partners and can be defined in terms of baryon, lepton and spin as $`R=(1)^{3B+L+2S}`$, is unbroken, there is at least one supersymmetric particle (the lightest supersymmetric particle or LSP) which must be stable. Thus, the minimal model contains the fewest number of new particles and interactions necessary to make a consistent theory.

There are very strong constraints, however, forbidding the existence of stable or long lived particles which are not color and electrically neutral . Strong and electromagnetically interacting LSPs would become bound with normal matter forming anomalously heavy isotopes. Indeed, there are very strong upper limits on the abundances, relative to hydrogen, of nuclear isotopes, $`n/n_H<10^{15}\mathrm{to}10^{29}`$ for 1 GeV $`<m<`$ 1 TeV. A strongly interacting stable relic is expected to have an abundance $`n/n_H<10^{10}`$ with a higher abundance for charged particles.

There are relatively few supersymmetric candidates which are not colored and are electrically neutral. The sneutrino is one possibility, but in the MSSM, it has been excluded as a dark matter candidate by direct and indirect searches. In fact, one can set an accelerator based limit on the sneutrino mass from neutrino counting, $`m_{\stackrel{~}{\nu }}>`$ 44.7 GeV . In this case, the direct relic searches in underground low-background experiments require $`m_{\stackrel{~}{\nu }}>`$ 20 TeV . Another possibility is the gravitino which is probably the most difficult to exclude. I will concentrate on the remaining possibility in the MSSM, namely the neutralinos.

## 2 Parameters

The most general version of the MSSM, despite its minimality in particles and interactions contains well over a hundred new parameters. The study of such a model would be untenable were it not for some (well motivated) assumptions. These have to do with the parameters associated with supersymmetry breaking. It is often assumed that, at some unification scale, all of the gaugino masses receive a common mass, $`m_{1/2}`$. The gaugino masses at the weak scale are determined by running a set of renormalization group equations. Similarly, one often assumes that all scalars receive a common mass, $`m_0`$, at the GUT scale. These too are run down to the weak scale. The remaining parameters of importance involve the Higgs sector. There is the Higgs mixing mass parameter, $`\mu `$, and since there are two Higgs doublets in the MSSM, there are two vacuum expectation values. One combination of these is related to the $`Z`$ mass, and therefore is not a free parameter, while the other combination, the ratio of the two vevs, $`\mathrm{tan}\beta `$, is free.

If the supersymmetry breaking Higgs soft masses are also unified at the GUT scale (and take the common value $`m_0`$), then $`\mu `$ and the physical Higgs masses at the weak scale are determined by electroweak vacuum conditions. This scenario is often referred to as the constrained MSSM or CMSSM. Once these parameters are set, the entire spectrum of sparticle masses at the weak scale can be calculated.

## 3 Neutralinos

There are four neutralinos, each of which is a linear combination of the $`R=1`$ neutral fermions,: the wino $`\stackrel{~}{W}^3`$, the partner of the 3rd component of the $`SU(2)_L`$ gauge boson; the bino, $`\stackrel{~}{B}`$, the partner of the $`U(1)_Y`$ gauge boson; and the two neutral Higgsinos, $`\stackrel{~}{H}_1`$ and $`\stackrel{~}{H}_2`$. Assuming gaugino mass universality at the GUT scale, the identity and mass of the LSP are determined by the gaugino mass $`m_{1/2}`$, $`\mu `$, and $`\mathrm{tan}\beta `$. In general, neutralinos can be expressed as a linear combination

$$\chi =\alpha \stackrel{~}{B}+\beta \stackrel{~}{W}^3+\gamma \stackrel{~}{H}_1+\delta \stackrel{~}{H}_2$$
(3)
The solution for the coefficients $`\alpha ,\beta ,\gamma `$ and $`\delta `$ for neutralinos that make up the LSP can be found by diagonalizing the mass matrix

$$(\stackrel{~}{W}^3,\stackrel{~}{B},\stackrel{~}{H}_{}^{0}{}_{1}{}^{},\stackrel{~}{H}_{}^{0}{}_{2}{}^{})\left(\begin{array}{cccc}M_2& 0& \frac{g_2v_1}{\sqrt{2}}& \frac{g_2v_2}{\sqrt{2}}\\ 0& M_1& \frac{g_1v_1}{\sqrt{2}}& \frac{g_1v_2}{\sqrt{2}}\\ \frac{g_2v_1}{\sqrt{2}}& \frac{g_1v_1}{\sqrt{2}}& 0& \mu \\ \frac{g_2v_2}{\sqrt{2}}& \frac{g_1v_2}{\sqrt{2}}& \mu & 0\end{array}\right)\left(\begin{array}{c}\stackrel{~}{W}^3\\ \stackrel{~}{B}\\ \stackrel{~}{H}_{}^{0}{}_{1}{}^{}\\ \stackrel{~}{H}_{}^{0}{}_{2}{}^{}\end{array}\right)$$
(4)
where $`M_1(M_2)`$ is a soft supersymmetry breaking term giving mass to the U(1) (SU(2)) gaugino(s). In a unified theory $`M_1=M_2`$ at the unification scale (at the weak scale, $`M_1\frac{5}{3}\frac{\alpha _1}{\alpha _2}M_2`$). As one can see, the coefficients $`\alpha ,\beta ,\gamma ,`$ and $`\delta `$ depend only on $`m_{1/2}`$, $`\mu `$, and $`\mathrm{tan}\beta `$.

In Figure 1 , regions in the $`M_2,\mu `$ plane with $`\mathrm{tan}\beta =2`$ are shown in which the LSP is one of several nearly pure states, the photino, $`\stackrel{~}{\gamma }`$, the bino, $`\stackrel{~}{B}`$, a symmetric combination of the Higgsinos, $`\stackrel{~}{H}_{(12)}`$, or the Higgsino, $`\stackrel{~}{S}=\mathrm{sin}\beta \stackrel{~}{H}_1+\mathrm{cos}\beta \stackrel{~}{H}_2`$. The dashed lines show the LSP mass contours. The cross hatched regions correspond to parameters giving a chargino ($`\stackrel{~}{W}^\pm ,\stackrel{~}{H}^\pm `$) state with mass $`m_{\stackrel{~}{\chi }}45GeV`$ and as such are excluded by LEP. This constraint has been extended by LEP and is shown by the light shaded region and corresponds to regions where the chargino mass is $`<104`$ GeV. The newer limit does not extend deep into the Higgsino region because of the degeneracy between the chargino and neutralino. Notice that the parameter space is dominated by the $`\stackrel{~}{B}`$ or $`\stackrel{~}{H}_{12}`$ pure states and that the photino only occupies a small fraction of the parameter space, as does the Higgsino combination $`\stackrel{~}{S}`$. Both of these light states are experimentally excluded.

## 4 The Relic Density

The relic abundance of LSP’s is determined by solving the Boltzmann equation for the LSP number density in an expanding Universe. The technique used is similar to that for computing the relic abundance of massive neutrinos. The relic density depends on additional parameters in the MSSM beyond $`m_{1/2},\mu `$, and $`\mathrm{tan}\beta `$. These include the sfermion masses, $`m_{\stackrel{~}{f}}`$ and the Higgs pseudo-scalar mass, $`m_A`$222In general, the relic density depends on the supersymmetry-breaking tri-linear masses $`A`$ (also assumed to be unified at the GUT scale) as well as two phases $`\theta _\mu `$ and $`\theta _A`$., derived from $`m_0`$ (and $`m_{1/2}`$). To determine the relic density it is necessary to obtain the general annihilation cross-section for neutralinos. In much of the parameter space of interest, the LSP is a bino and the annihilation proceeds mainly through sfermion exchange. Because of the p-wave suppression associated with Majorana fermions, the s-wave part of the annihilation cross-section is suppressed by the outgoing fermion masses. This means that it is necessary to expand the cross-section to include p-wave corrections which can be expressed as a term proportional to the temperature if neutralinos are in equilibrium. Unless the neutralino mass happens to lie near near a pole, such as $`m_\chi `$ $`m_Z/2`$ or $`m_h/2`$, in which case there are large contributions to the annihilation through direct $`s`$-channel resonance exchange, the dominant contribution to the $`\stackrel{~}{B}\stackrel{~}{B}`$ annihilation cross section comes from crossed $`t`$-channel sfermion exchange.

Annihilations in the early Universe continue until the annihilation rate $`\mathrm{\Gamma }\sigma vn_\chi `$ drops below the expansion rate given by the Hubble parameter, $`H`$. For particles which annihilate through approximate weak scale interactions, this occurs when $`Tm_\chi /20`$. Subsequently, the relic density of neutralinos is fixed relative to the number of relativistic particles. As noted above, the number density of neutralinos is tracked by a Boltzmann-like equation,

$$\frac{dn}{dt}=3\frac{\dot{R}}{R}n\sigma v(n^2n_0^2)$$
(5)
where $`n_0`$ is the equilibrium number density of neutralinos. By defining the quantity $`f=n/T^3`$, we can rewrite this equation in terms of $`x`$, as

$$\frac{df}{dx}=m_\chi \left(\frac{1}{90}\pi ^2\kappa ^2N\right)^{1/2}(f^2f_0^2)$$
(6)
The solution to this equation at late times (small $`x`$) yields a constant value of $`f`$, so that $`nT^3`$. The final relic density expressed as a fraction of the critical energy density can be written as

$$\mathrm{\Omega }_\chi h^21.9\times 10^{11}\left(\frac{T_\chi }{T_\gamma }\right)^3N_f^{1/2}\left(\frac{\mathrm{GeV}}{ax_f+\frac{1}{2}bx_f^2}\right)$$
(7)
where $`(T_\chi /T_\gamma )^3`$ accounts for the subsequent reheating of the photon temperature with respect to $`\chi `$, due to the annihilations of particles with mass $`m<x_fm_\chi `$ . The subscript $`f`$ refers to values at freeze-out, i.e., when annihilations cease. The coefficients $`a`$ and $`b`$ are related to the partial wave expansion of the cross-section, $`\sigma v=a+bx+\mathrm{}`$. Eq. (7 ) results in a very good approximation to the relic density expect near s-channel annihilation poles, thresholds and in regions where the LSP is nearly degenerate with the next lightest supersymmetric particle.

## 5 Phenomenological and Cosmological Constraints

For the cosmological limits on the relic density I will assume

$$0.1\mathrm{\Omega }_\chi h^2\mathrm{\hspace{0.33em}0.3}.$$
(8)
The upper limit being a conservative bound based only on the lower limit to the age of the Universe of 12 Gyr. Indeed, most analyses indicate that $`\mathrm{\Omega }_{\mathrm{matter}}<0.40.5`$ and thus it is very likely that $`\mathrm{\Omega }_\chi h^2<0.2`$. One should note that smaller values of $`\mathrm{\Omega }_\chi h^2`$ are allowed, since it is quite possible that some of the cold dark matter might not consist of LSPs.

The calculated relic density is found to have a relevant cosmological density over a wide range of susy parameters. For all values of $`\mathrm{tan}\beta `$, there is a ‘bulk’ region with relatively low values of $`m_{1/2}`$ and $`m_0`$ where $`0.1<\mathrm{\Omega }_\chi h^2<0.3`$. However there are a number of regions at large values of $`m_{1/2}`$ and/or $`m_0`$ where the relic density is still compatible with the cosmological constraints. At large values of $`m_{1/2}`$, the lighter stau, becomes nearly degenerate with the neutralino and co-annihilations between these particles must be taken into account. For non-zero values of $`A_0`$, there are new regions for which $`\chi \stackrel{~}{t}`$ coannihilations are important. At large $`\mathrm{tan}\beta `$, as one increases $`m_{1/2}`$, the pseudo-scalar mass, $`m_A`$ begins to drop so that there is a wide funnel-like region (at all values of $`m_0`$) such that $`2m_\chi m_A`$ and s-channel annihilations become important. Finally, there is a region at very high $`m_0`$ where the value of $`\mu `$ begins to fall and the LSP becomes more Higgsino-like. This is known as the ‘focus point’ region.

As an aid to the assessment of the prospects for detecting sparticles at different accelerators, benchmark sets of supersymmetric parameters have often been found useful, since they provide a focus for concentrated discussion. A set of proposed post-LEP benchmark scenarios in the CMSSM are illustrated schematically in Fig. 2. Five of the chosen points are in the ‘bulk’ region at small $`m_{1/2}`$ and $`m_0`$, four are spread along the coannihilation ‘tail’ at larger $`m_{1/2}`$ for various values of $`\mathrm{tan}\beta `$. This tail runs along the shaded region in the lower right corner where the stau is the LSP and is therefore excluded by the constraints against charged dark matter. Two points are in rapid-annihilation ‘funnels’ at large $`m_{1/2}`$ and $`m_0`$. At large values of $`m_0`$, the focus-point region runs along the boundary where electroweak symmetry no longer occurs (shown in Fig. 2 as the shaded region in the upper left corner). Two points were chosen in the focus-point region at large $`m_0`$. The proposed points range over the allowed values of $`\mathrm{tan}\beta `$ between 5 and 50. The light shaded region corresponds to the portion of parameter space where the relic density $`\mathrm{\Omega }_\chi h^2`$ is between 0.1 and 0.3.

The most important phenomenological constraints are also shown schematically in Figure 2. These include the constraint provided by the LEP lower limit on the Higgs mass: $`m_H>`$ 114.1 GeV . This holds in the Standard Model, for the lightest Higgs boson $`h`$ in the general MSSM for $`\mathrm{tan}\beta <8`$, and almost always in the CMSSM for all $`\mathrm{tan}\beta `$. Since $`m_h`$ is sensitive to sparticle masses, particularly $`m_{\stackrel{~}{t}}`$, via loop corrections, the Higgs limit also imposes important constraints on the CMSSM parameters, principally $`m_{1/2}`$ as seen by the dashed curve in Fig. 2. The constraint imposed by measurements of $`bs\gamma `$  also exclude small values of $`m_{1/2}`$. These measurements agree with the Standard Model, and therefore provide bounds on MSSM particles, such as the chargino and charged Higgs masses, in particular. Typically, the $`bs\gamma `$ constraint is more important for $`\mu <0`$, but it is also relevant for $`\mu >0`$, particularly when $`\mathrm{tan}\beta `$ is large. The BNL E821 experiment reported last year a new measurement of $`a_\mu \frac{1}{2}(g_\mu 2)`$ which deviated by 2.6 standard deviations from the best Standard Model prediction available at that time . The largest contribution to the errors in the comparison with theory was thought to be the statistical error of the experiment, which has been significantly reduced just recently . However, it has recently been realized that the sign of the most important pseudoscalar-meson pole part of the light-by-light scattering contribution to the Standard Model prediction should be reversed, which reduces the apparent experimental discrepancy to about 1.6 standard deviations ($`\delta a_\mu \times 10^{10}=26\pm 16`$). With the new data, the discrepancy with theory ranges from 1.6 to 2.6 $`\sigma `$, i.e., $`\delta a_\mu \times 10^{10}=26\pm 10`$ to $`17\pm 11`$. This constraint excludes very small values of $`m_{1/2}`$ and $`m_0`$. In Fig. 2, the $`g2`$ constraint is shown schematically by the dotted line. It may also exclude very large values of the parameters as well as negative values of $`\mu `$, if the discrepancy holds up.

Following a previous analysis, in Figure 3 the $`m_{1/2}m_0`$ parameter space is shown for $`\mathrm{tan}\beta =10`$. The dark shaded region (in the lower right) corresponds to the parameters where the LSP is not a neutralino but rather a $`\stackrel{~}{\tau }_R`$. The cosmologically interesting region at the left of the figure is due to the appearance of pole effects. There, the LSP can annihilate through s-channel $`Z`$ and $`h`$ (the light Higgs) exchange, thereby allowing a very large value of $`m_0`$. However, this region is excluded by phenomenological constraints. Here one can see clearly the coannihilation tail which extends towards large values of $`m_{1/2}`$. In addition to the phenomenological constraints discussed above, Figure 3 also shows the current experimental constraints on the CMSSM parameter space due to the limit $`m_{\chi ^\pm }>`$ 103.5 GeV provided by chargino searches at LEP . LEP has also provided lower limits on slepton masses, of which the strongest is $`m_{\stackrel{~}{e}}>`$ 99 GeV . This is shown by dot-dashed curve in the lower left corner of Fig. 3.

As one can see, one of the most important phenomenological constraint at this value of $`\mathrm{tan}\beta `$ is due to the Higgs mass (shown by the nearly vertical dot-dashed curve). The theoretical Higgs masses were evaluated using FeynHiggs, which is estimated to have a residual uncertainty of a couple of GeV in $`m_h`$. The region excluded by the $`bs\gamma `$ constraint is the dark shaded (green) region to the left of the plot.

As many authors have pointed out, a discrepancy between theory and the BNL experiment could well be explained by supersymmetry. As seen in Fig. 3, this is particularly easy if $`\mu >0`$. The medium (pink) shaded region in the figure corresponds to the overall allowed region by the new experimental result: $`5<\delta a_\mu \times 10^{10}<46`$. The two solid lines within the shaded region corresponds to the central values $`\delta a_\mu \times 10^{10}=17`$ and 26 respectively. The optimistic $`2\sigma `$ lower bound of $`\delta a_\mu \times 10^{10}=6`$ is shown as a dashed curve.

As discussed above, another mechanism for extending the allowed CMSSM region to large $`m_\chi `$ is rapid annihilation via a direct-channel pole when $`m_\chi \frac{1}{2}m_A`$ . This may yield a ‘funnel’ extending to large $`m_{1/2}`$ and $`m_0`$ at large $`\mathrm{tan}\beta `$, as seen in Fig. 4.

In principle the true input parameters in the CMSSM are: $`\mu ,m_1,m_2,`$ and $`B`$, where $`m_1`$ and $`m_2`$ are the Higgs soft masses (in the CMSSM $`m_1=m_2=m_0`$ and $`B`$ is the susy breaking bilinear mass term). In this case, the electroweak symmetry breaking conditions lead to a prediction of $`M_Z,\mathrm{tan}\beta `$ ,and $`m_A`$. Since we are not really interested in predicting $`M_Z`$, it is more useful to assume instead the following CMSSM input parameters: $`M_Z,m_1,m_2,`$ and $`\mathrm{tan}\beta `$ again with $`m_1=m_2=m_0`$. In this case, one predicts $`\mu ,B`$, and $`m_A`$. However, one can generalize the CMSSM case to include non-universal Higgs masses (NUHM), in which case the input parameters become:$`M_Z,\mu ,m_A,`$ and $`\mathrm{tan}\beta `$ and one predicts $`m_1,m_2`$, and $`B`$.

The NUHM parameter space was recently analyzed and a sample of the results found is shown in Fig. 5. While much of the cosmologically preferred area with $`\mu <0`$ is excluded, there is a significant enhancement in the allowed parameter space for $`\mu >0`$.

### 5.1 Detection

Because the LSP as dark matter is present locally, there are many avenues for pursuing dark matter detection. Direct detection techniques rely on an ample neutralino-nucleon scattering cross-section. The effective four-fermion lagrangian can be written as

$``$ $`=`$ $`\overline{\chi }\gamma ^\mu \gamma ^5\chi \overline{q_i}\gamma _\mu (\alpha _{1i}+\alpha _{2i}\gamma ^5)q_i`$ (9)
$`+`$ $`\alpha _{3i}\overline{\chi }\chi \overline{q_i}q_i+\alpha _{4i}\overline{\chi }\gamma ^5\chi \overline{q_i}\gamma ^5q_i`$
$`+`$ $`\alpha _{5i}\overline{\chi }\chi \overline{q_i}\gamma ^5q_i+\alpha _{6i}\overline{\chi }\gamma ^5\chi \overline{q_i}q_i`$
However, the terms involving $`\alpha _{1i},\alpha _{4i},\alpha _{5i}`$, and $`\alpha _{6i}`$ lead to velocity dependent elastic cross sections. The remaining terms are: the spin dependent coefficient, $`\alpha _{2i}`$ and the scalar coefficient $`\alpha _{3i}`$. Contributions to $`\alpha _{2i}`$ are predominantly through light squark exchange. This is the dominant channel for binos. Scattering also occurs through Z exchange but this channel requires a strong Higgsino component. Contributions to $`\alpha _{3i}`$ are also dominated by light squark exchange but Higgs exchange is non-negligible in most cases.

The results from a CMSSM and MSSM analysis for $`\mathrm{tan}\beta =3`$ and 10 are compared with the most recent CDMS and Edelweiss bounds in Fig. 6. These results have nearly entirely excluded the region purported by the DAMA experiment. The CMSSM prediction  is shown by the dark shaded region, while the NUHM case  is shown by the larger lighter shaded region.

I conclude by showing the prospects for direct detection for the benchmark points discussed above. Fig. 7 shows rates for the elastic spin-independent scattering of supersymmetric relics, including the projected sensitivities for CDMS II  and CRESST (solid) and GENIUS (dashed). Also shown are the cross sections calculated in the proposed benchmark scenarios discussed in the previous section, which are considerably below the DAMA range ($`10^510^6`$ pb). Indirect searches for supersymmetric dark matter via the products of annihilations in the galactic halo or inside the Sun also have prospects in some of the benchmark scenarios.

## Acknowledgments

I would like to thank J. Ellis, T. Falk, A. Ferstl, G. Ganis, Y. Santoso, and M. Srednicki for enjoyable collaborations from which this work is culled. This work was supported in part by DOE grant DE-FG02-94ER40823 at Minnesota.


text:
# QED Radiative Corrections in Processes of Exclusive Pion Electroproduction

## I Introduction

Understanding the electromagnetic transition amplitudes from the nucleon ground state to excited states provides valuable insights into the electromagnetic structure of the nucleon. Exclusive pion electroproduction is one of the major sources to provide the most direct information about the spatial and spin structure of the excited states. With the development of high intensity and high duty–factor electron beam with high degree of polarization, this field reaches a new level of quality. For the past several years, exclusive pion electroproduction has been the main subject of extensive studies at various accelerator laboratories such as MIT-Bates, ELSA, MAMI and Jefferson Lab.

New measurements with the CLAS detector at Jefferson Lab/Hall B nstar  are expected to greatly improve the systematic and statistical precision and cover a wide kinematic range in four–momentum transfer $`Q^2`$ and invariant mass $`W`$, as well as the full angular range of the resonance decay into the nucleon–pion final state:

$`e(k_1)+p(p)`$ $``$ $`e^{}(k_2)+\pi ^+(p_h)+n(p_u),`$
$`e(k_1)+p(p)`$ $``$ $`e^{}(k_2)+p(p_h)+\pi ^0(p_u),`$ (1)
where $`p_h(p_u)`$ denotes the momentum of the detected (undetected) hadron. Adequate calculation of radiative corrections (RC) becomes important in interpreting the measured observables such as unpolarized coincidence cross sections and polarization asymmetries.

While solving the RC problem, people are commonly referred to the classical approach developed by Mo and Tsai MoTsai  and used for inclusive and elastic electron scattering for decades. However, this approach cannot be directly applied for exclusive pion electroproduction due to the following reasons.

First, we now deal with exclusive electroproduction, where the hadron is detected in addition to the final electron. It reduces the room for phase space allowed for the final radiated photon. The formulas of Mo and Tsai as well as any other inclusive formulas cannot be applied for this case without additional strong assumptions.

Second, there are only contributions of two structure functions in the inclusive case. The exact formalism for the exclusive process requires consideration of four structure functions for the unpolarized case with additional angular dependence associated with them. Note that transition to the case of two structure functions would not be possible even within realistic approximations. Mo and Tsai’s approach predicts neither RC to polarization asymmetries, nor dependencies on the outgoing hadron angles.

The third reason is a known shortcoming of the Mo and Tsai’s approach, namely, the dependence on an unphysical parameter splitting soft and hard regions of the phase space of radiated photon in order to cancel the infrared divergence.

In our approach, which is based on a covariant procedure of infrared divergence cancellation proposed by Bardin and Shumeiko in Ref.BSh , such an unphysical parameter is not required. Previously, this approach was applied for the calculation of RC for inclusive ABSh ; KSh1 ; ASh , semi–inclusive Soroko ; AkuSor  and exclusive diffractive Ak  reactions. Recent reviews of the approach, higher order effects and calculation for specific experiments can be found in papers ABK ; AKS ; AISh2001 . Based on these results, a Fortran code POLRAD POLRAD20  for RC calculation in polarized inclusive and semi–inclusive processes was developed. Besides, some specific tasks such as Monte Carlo generator RADGEN RADGEN , MC approach to diffractive vector meson electroproduction DIFF2 , RC to spin–density matrix elements in exclusive vector meson production AKuzh , quasielastic tail for polarized He–3 target ARST  were solved. Recently Bardin–Shumeiko approach was applied to the measurements of elastic polarized electron–proton scattering at Jefferson Lab pep1 ; pep2 . A comprehensive analysis of results obtained in Ref. MoTsai  and BSh  was made in Refs. RADGEN ; KuzhSh .

Note that there are also other approaches for calculation of RC in electroproduction processes. For example, the results of Ref. spiesberger  are actively used for HERA experiments, while the approach developed in papers Kuraev  is applied to specific measurements at Jefferson Lab AAMsf .

The Feynman diagrams needed to calculate RC are presented in Fig.1. They include QED processes of radiation of an unobserved real photon, vacuum polarization and lepton–photon vertex corrections. These processes give the largest contribution due to a large logarithmic term $`\mathrm{ln}(Q^2/m_e^2)`$. They can be calculated exactly from QED, and uncertainties of such a calculation are only due to the fits and data used for the hadronic structure functions. These uncertainties are demonstrated in the present article. Additional mechanisms (box–type diagrams, emission by hadrons) are smaller by about an order of magnitude and they contain considerable theoretical uncertainties. Most recent studies of two–photon exchange effects for elastic $`ep`$–scattering were reported in Refs.MaxTj ; AAExcl02 . Previous RC calculations AISh95  for inclusive polarized deep–inelastic scattering included above–mentioned additional mechanisms. Generalization to the case of meson electroproduction will be subject to a separate study.

The paper is organized as follows. We introduce kinematics and definitions (Section II), derive a cross section of the radiative process (Section III), solve the problem of infrared divergence (Section IV), obtain RC in the leading log approximation (Section V), verify the relation between exclusive and inclusive RC (Section VI), perform numerical analysis (Section VII), and summarize the results in Section VIII.

## II Kinematics and Born process

At the Born level (Fig.1a), the cross section of the processes (I) is described by four kinematic variables. Following tradition, we choose them as squared virtual photon momentum $`Q^2`$, invariant mass of initial proton and the virtual photon $`W`$, and detected pion (or proton) angles $`\theta _h`$ and $`\varphi _h`$ in the center–of–mass of the final hadrons. The obtained formulas are equally valid for both electron and muon scattering.

We use the following Lorentz invariants defined from leptonic 4–momenta:

$`S=2k_1p,X=2k_2p,Q^2=(k_1k_2)^2,`$
$`u_1=SQ^2,u_2=X+Q^2,\lambda _{1,2}=u_{1,2}^24m^2W^2,`$
$`S_{p,x}=S\pm X,\lambda _S=S^24m^2M^2`$
$`W^2=S_xQ^2+M^2,\lambda _q=S_x^2+4M^2Q^2,`$
$`\lambda =Q^2u_1u_2Q^4W^2m^2\lambda _q,`$ (2)
where $`m(M)`$ is the lepton (proton) mass.

We use the c.m. system of virtual photon ($`q=k_1k_2`$) and initial nucleon. The axis $`OZ`$ is chosen along $`\stackrel{}{q}`$ (see Fig. 2). The energies and angles (Fig.2) in the selected frame can be expressed in terms of invariants (II) as follows:

$`E_{1,2}={\displaystyle \frac{u_{1,2}}{2W}},E_q={\displaystyle \frac{S_x2Q^2}{2W}},E_p={\displaystyle \frac{S_x+2M^2}{2W}},`$
$`p_{1,2}={\displaystyle \frac{\sqrt{\lambda _{1,2}}}{2W}},p_p=p_q={\displaystyle \frac{\sqrt{\lambda _q}}{2W}},`$
$`\mathrm{cos}\theta _{1,2}={\displaystyle \frac{u_{1,2}(S_x2Q^2)\pm 2Q^2W^2}{\sqrt{\lambda _{1,2}}\sqrt{\lambda _q}}},`$
$`\mathrm{sin}\theta _{1,2}={\displaystyle \frac{2W\sqrt{\lambda }}{\sqrt{\lambda _{1,2}}\sqrt{\lambda _q}}}.`$ (3)
Here subscripts $`1,2,q,p`$ denote the initial lepton, final lepton, virtual photon and initial proton, respectively.

All the introduced kinematic variables have the same definition for both Born and radiative kinematics. But the situation is different for the final hadrons. The reason is that the energy of the observed hadron is not fixed by measurements of the chosen kinematic variables: $`Q^2`$, $`W^2`$, $`\theta _h`$ $`\varphi _h`$. As a result, it is different for these cases. In the Born case the hadron energy can be defined via conservation laws, but for the radiative process it depends on the unknown energy of the radiated photon. In order to distinguish between the variables, we use the superscript ’0’ for the quantities calculated for Born kinematics:

$`E_h^0={\displaystyle \frac{W^2+m_h^2m_u^2}{2W}},p_h^0={\displaystyle \frac{\sqrt{\lambda _W^0}}{2W}},`$
$`\lambda _W^0=(W^2m_h^2m_u^2)^24m_h^2W^2.`$ (4)
All invariants defined via the measured 4–momentum of hadron $`p_h`$ are also different for Born and radiative cases:

$`V_{1,2}^0`$ $`=`$ $`2k_{1,2}p_h=2(E_{1,2}E_h^0p_{1,2}p_h^0(\mathrm{cos}\theta _h\mathrm{cos}\theta _{1,2}+`$
$`\mathrm{sin}\theta _h\mathrm{sin}\theta _{1,2}\mathrm{cos}\varphi _h)),`$
$`S_t^0`$ $`=`$ $`2pp_h=W^2+m_h^2m_u^2+V_2^0V_1^0.`$ (5)
The last expression follows from 4–momentum conservation.

The Born cross section of exclusive pion electroproduction in terms of the introduced variables reads:

$`d\sigma _0`$ $`=`$ $`{\displaystyle \frac{M_0^2}{2S(2\pi )^5}}{\displaystyle \frac{d\stackrel{}{k}_2}{2E_2}}{\displaystyle \frac{d\stackrel{}{p}_u}{2E_u}}{\displaystyle \frac{d\stackrel{}{p}_h}{2E_h}}\delta (\mathrm{\Lambda }p_u)=`$ (6)
$`={\displaystyle \frac{(4\pi \alpha )^2}{2(4\pi )^4S^2}}{\displaystyle \frac{L_{\mu \nu }^0W_{\mu \nu }}{Q^4}}{\displaystyle \frac{\sqrt{\lambda _W^0}}{W^2}}dQ^2dW^2d\mathrm{\Omega }_h,`$
where $`\mathrm{\Lambda }`$ is the total 4–momentum of the undetected particles, $`\mathrm{\Lambda }=p+k_1k_2p_h`$. For the Born process (Fig.1a), $`\mathrm{\Lambda }`$ is equal to the momentum of the undetected hadron, $`p_u`$, while for the radiative process (Fig.1b,c), it is a sum of the undetected hadron and bremsstrahlung photon momenta, $`\mathrm{\Lambda }=p_u+k`$.

The phase space is calculated as

$$\frac{d\stackrel{}{k}_2}{2E_2}=\frac{\pi }{2S}dQ^2dW^2$$
(7)
and

$$\frac{d\stackrel{}{p}_u}{2E_u}\frac{d\stackrel{}{p}_h}{2E_h}\delta (\mathrm{\Lambda }p_u)=\frac{p_h^0}{4W}d\mathrm{\Omega }_h=\frac{\sqrt{\lambda _W^0}}{8W^2}d\mathrm{\Omega }_h.$$
(8)
To calculate the matrix element squared, one needs to contract the leptonic and hadronic tensors,

$$M_0^2=\frac{e^4}{Q^4}L_{\mu \nu }^0W_{\mu \nu }=\frac{2e^4}{Q^4}\underset{i=1}{\overset{5}{}}\theta _i^0_i^0.$$
(9)
We consider the longitudinally polarized lepton beam. In this case the leptonic tensor reads:

$$L_{\mu \nu }^0=\frac{1}{2}\mathrm{Tr}(\widehat{k}_2+m)\gamma _\mu (\widehat{k}_1+m)(1+i\gamma _5\widehat{\xi })\gamma _\nu .$$
(10)
Here the lepton polarization vector is kept in a general form. If the lepton is longitudinally polarized and its helicity is positive, then the vector $`\xi `$ can be expressed as ASh

$$\xi =\frac{1}{\sqrt{\lambda _s}}\left(\frac{S}{m}k_12mp\right).$$
(11)
In the Born approximation, the second term can be dropped. However, it gives a non–zero contribution to RC.

For the hadronic tensor, we use a general covariant form

$`W^{\mu \nu }`$ $`=`$ $`\stackrel{~}{g}^{\mu \nu }_1+\stackrel{~}{p}^\mu \stackrel{~}{p}^\nu _2+\stackrel{~}{p}_{h}^{}{}_{}{}^{\mu }\stackrel{~}{p}_{h}^{}{}_{}{}^{\nu }_3`$
$`+(\stackrel{~}{p}^\mu \stackrel{~}{p}_{h}^{}{}_{}{}^{\nu }+\stackrel{~}{p}_{h}^{}{}_{}{}^{\mu }\stackrel{~}{p}^\nu )_4+(\stackrel{~}{p}_{h}^{}{}_{}{}^{\mu }\stackrel{~}{p}^\nu \stackrel{~}{p}^\mu \stackrel{~}{p}_{h}^{}{}_{}{}^{\nu })_5,`$
where the tilde for an arbitrary 4–vector $`a^\mu `$ denotes the substitution $`\stackrel{~}{a}^\mu =a^\mu \frac{aq}{q^2}q^\mu `$ (to ensure electromagnetic gauge invariance). The first four structure functions have tensor coefficients symmetric over Lorentz indices, but the last one is antisymmetric, it contributes to the polarization–dependent part of the cross section. The quantities $`\theta _i^0`$ have the following form,

$`\theta _1^0`$ $`=`$ $`Q^2,`$ (13)
$`\theta _2^0`$ $`=`$ $`{\displaystyle \frac{1}{2}}(SXM^2Q^2),`$
$`\theta _3^0`$ $`=`$ $`{\displaystyle \frac{1}{2}}(V_1^0V_2^0m_h^2Q^2),`$
$`\theta _4^0`$ $`=`$ $`{\displaystyle \frac{1}{2}}(SV_2^0+XV_1^0S_t^0Q^2),`$
$`\theta _5^0`$ $`=`$ $`2ϵ(k_1,k_2,p,p_h),`$
where the operator in the last line is defined as

$$ϵ(p_1,p_2,p_3,p_4)=p_1^\alpha p_2^\beta p_3^\gamma p_4^\sigma ϵ_{\alpha \beta \gamma \sigma }.$$
(14)
In the chosen c.m. frame it is equal to

$`ϵ(k_1,k_2,p,p_h)`$ $`=`$ $`WP_qp_hp_1\mathrm{sin}\theta _1\mathrm{sin}\theta _h\mathrm{sin}\varphi _h`$ (15)
$`=`$ $`{\displaystyle \frac{\sqrt{\lambda _W^0\lambda }}{4W}}\mathrm{sin}\theta _h\mathrm{sin}\varphi _h.`$
As a result, we obtain for the Born cross section, ($`\sigma _0=d\sigma /dW^2dQ^2d\mathrm{\Omega }_h`$)

$$\sigma _0=\frac{\alpha ^2\sqrt{\lambda _W^0}}{32\pi ^2S^2W^2Q^4}\underset{i=1}{\overset{5}{}}\theta _i^0_i^0.$$
(16)
The expression for the Born cross section in the form (16) is convenient for our further calculation. It is equivalent to the well known formula in terms of photoabsorption cross sections (or response functions), where each term corresponds to certain polarization states of the virtual photon dombey ,

$`{\displaystyle \frac{1}{N_0}}\sigma _0`$ $`=`$ $`\sigma _T+ϵ\sigma _L+ϵ\mathrm{cos}2\varphi _h\sigma _{TT}`$ (17)
$`+\sqrt{ϵ(1+ϵ)/2}\mathrm{cos}\varphi _h\sigma _{LT}`$
$`+h_e\sqrt{ϵ(1ϵ)/2}\mathrm{sin}\varphi _h\sigma _{LT}^{},`$
where

$$N_0=\alpha \frac{W^2M^2}{2S^2Q^2(1ϵ)},$$
(18)
$`h_e`$ is the longitudinal polarization degree of the incoming leptons, and $`ϵ`$ dombey  describes the virtual photon polarization.

The following expressions relate the structure functions (II) and Born photoabsorption cross sections (17),

$`_1(W^2,Q^2,t)`$ $`=`$ $`C(\sigma _T\sigma _{TT}),`$
$`_2(W^2,Q^2,t)`$ $`=`$ $`{\displaystyle \frac{2C}{\lambda _q}}(2Q^2(\sigma _T\sigma _{TT}+\sigma _L),`$
$`TQ\sigma _{LT}+T^2\sigma _{TT}),`$
$`_3(W^2,Q^2,t)`$ $`=`$ $`{\displaystyle \frac{2C\lambda _q}{\lambda _l}}\sigma _{TT},`$
$`_4(W^2,Q^2,t)`$ $`=`$ $`{\displaystyle \frac{C}{\sqrt{\lambda _l}}}(2T\sigma _{TT}Q\sigma _{LT}),`$
$`_5(W^2,Q^2,t)`$ $`=`$ $`{\displaystyle \frac{C}{\sqrt{\lambda _l}}}Q\sigma _{LT}^{},`$ (19)
where

$`C`$ $`=`$ $`{\displaystyle \frac{16\pi ^2(W^2M^2)W^2}{\alpha \sqrt{\lambda _W^0}}},`$
$`T`$ $`=`$ $`{\displaystyle \frac{S_xt_q2Q^2S_t}{\sqrt{\lambda _l}}},`$
$`\lambda _l`$ $`=`$ $`Q^2S_t^2S_tS_xt_qM^2t_q^2m_h^2\lambda _q,`$ (20)
and $`Q=\sqrt{Q^2}`$. Note that both structure functions $`_i`$, as well as cross sections $`\sigma `$’s, are functions of three independent invariant variables, which usually are chosen as $`Q^2`$, $`W^2`$ and $`t`$. Therefore their transformation coefficients depend only on these variables but not on $`\varphi _h`$. For the radiative case the variables will be defined in the next section (see Eq.(30)). For the Born case, they are taken as: $`W^2`$, $`Q^2`$ and $`t=t_0=V_2^0V_1^0Q^2+m_h^2`$.

## III Exact formulas for radiative correction

The cross section of radiative process is given by

$`d\sigma _r`$ $`=`$ $`{\displaystyle \frac{M_r^2}{2S(2\pi )^8}}{\displaystyle \frac{d\stackrel{}{k}_2}{2E_2}}{\displaystyle \frac{d\stackrel{}{k}}{2\omega }}{\displaystyle \frac{d\stackrel{}{p}_u}{2E_u}}{\displaystyle \frac{d\stackrel{}{p}_h}{2E_h}}\delta (\mathrm{\Lambda }kp_u)`$
$`=`$ $`{\displaystyle \frac{(4\pi \alpha )^3dQ^2dW^2d\mathrm{\Omega }_h}{2(4\pi )^7S^2W^2}}{\displaystyle 𝑑\mathrm{\Omega }_k𝑑v\frac{v\sqrt{\lambda _W}}{f_W^2\stackrel{~}{Q}^4}L_{\mu \nu }^RW_{\mu \nu }},`$
where $`\mathrm{\Omega }_{k(h)}`$ stands for the solid angle of the bremsstrahlung photon (detected hadron).

Let us first consider the phase space for the radiative process. Integrating over the 3-momentum of unobserved hadron and using Dirac $`\delta `$-function to eliminate integration over the hadron energy, we have

$`{\displaystyle \frac{d\stackrel{}{k}}{2\omega }\frac{d\stackrel{}{p}_u}{2E_u}\frac{d\stackrel{}{p}_h}{2E_h}\delta (\mathrm{\Lambda }kp_u)}=`$
$`={\displaystyle \frac{1}{64}}{\displaystyle 𝑑\mathrm{\Omega }_h𝑑\mathrm{\Omega }_k\underset{0}{\overset{v_m}{}}\frac{vdv}{f_W^2}\sqrt{\lambda _W}},`$ (22)
where

$`f_W`$ $`=`$ $`WE_h+p_h(\mathrm{cos}\theta _h\mathrm{cos}\theta _k`$ (23)
$`+\mathrm{sin}\theta _h\mathrm{sin}\theta _k\mathrm{cos}(\varphi _h\varphi _k)).`$
Here we introduce a quantity $`v`$ that describes the missing mass (or inelasticity) due to emission of a bremsstrahlung photon, $`v=\mathrm{\Lambda }^2m_u^2`$. Note that $`v=0`$ for the Born process (Fig.1a), as well as in the soft–photon limit ($`k=0`$). As can be seen from (III), both the photon and hadron energies are now functions of $`v`$. It is related to another quantity, $`f_W`$, as follows:

$$v=W^2+m_h^2m_u^22WE_h=2\omega f_W.$$
(24)
The largest value of inelasticity allowed by kinematics ($`v_m`$) corresponds to the threshold of electroproduction. It is therefore defined from the relation $`E_h=m_h`$, yielding

$$v_m=(Wm_h)^2m_u^2.$$
(25)
Note that $`v_m`$ is always smaller for the heavier hardon, namely, the nucleon, detected in the final state. It does not depend on photon angles, therefore the integration region is a rectangle.

The maximum inelasticity $`v_m`$ is an important quantity for the RC calculation. All kinematic cuts made by experimentalists in data analysis influence RC. It is often possible to reduce all these cuts to one effective cut on the inelasticity (or missing mass), $`v_{cut}`$. In which case $`v_{cut}`$ should replace $`v_m`$ as the upper limit of integration in Eq.(III), and thus RC can be calculated within the cuts using the obtained formulas. If no cuts are applied, the maximum value of $`v`$ equals $`v_m`$, as given by energy–momentum conservation (25).

Now we can fix the kinematics of the radiative process. We have to express all scalar products and kinematic variables in terms of seven variables: Four variables that define the differential cross section and three integration variables. First, note that all definitions of leptonic variables given in (II) hold in this case. Hadron and real photon energies are defined by (24). Hadron momentum and scalar products are given as

$`p_h`$ $`=`$ $`{\displaystyle \frac{\sqrt{\lambda _W}}{2W}},\lambda _W=(W^2m_h^2m_u^2v)^24m_h^2W^2,`$
$`V_{1,2}`$ $`=`$ $`2k_{1,2}p_h=2(E_{1,2}E_hp_{1,2}p_h(\mathrm{cos}\theta _h\mathrm{cos}\theta _{1,2}+`$
$`\mathrm{sin}\theta _h\mathrm{sin}\theta _{1,2}\mathrm{cos}\varphi _h)),`$
$`S_t`$ $`=`$ $`2pp_h=W^2+m_h^2m_u^2+V_2V_1v`$
$`=S_x+t+M^2m_u^2v,`$
$`t`$ $`=`$ $`V_2V_1Q^2+m_h^2.`$ (26)
As in (II), the last expression is obtained from conservation laws. All hadron kinematic variables depend on inelasticity $`v`$ but not on photon angles $`\mathrm{\Omega }_k`$. Scalar products containing the photon 4–momentum read

$`2pk`$ $`=`$ $`R_w(1\tau )=2\omega (E_pp_p\mathrm{cos}\theta _k),`$ (27)
$`2k_{1,2}k`$ $`=`$ $`\kappa _{1,2}=2\omega (E_{1,2}p_{1,2}(\mathrm{cos}\theta _k\mathrm{cos}\theta _{1,2}+`$
$`\mathrm{sin}\theta _k\mathrm{sin}\theta _{1,2}\mathrm{cos}\varphi _k),`$
$`2p_hk`$ $`=`$ $`\mu =2\omega (E_hp_h(\mathrm{cos}\theta _k\mathrm{cos}\theta _h+`$
$`\mathrm{sin}\theta _k\mathrm{sin}\theta _h\mathrm{cos}(\varphi _k\varphi _h)).`$
Here we use also invariant variables $`R_w=2k(p+q)`$ and $`\tau =2kq/R_w`$; $`R_w=f_Wv/W`$.

The leptonic tensor of the radiative process has a more complicated form. It can be written as

$`L_{\mu \nu }^R`$ $`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{Tr}(\widehat{k}_2+m)\mathrm{\Gamma }_{\mu \alpha }(\widehat{k}_1+m)(1+i\gamma _5\widehat{\xi })\widehat{\mathrm{\Gamma }}_{\alpha \nu },`$
$`\mathrm{\Gamma }_{\mu \alpha }`$ $`=`$ $`\left[\left({\displaystyle \frac{k_{1\alpha }}{kk_1}}{\displaystyle \frac{k_{1\alpha }}{kk_1}}\right)\gamma _\mu {\displaystyle \frac{\gamma _\mu \widehat{k}\gamma \alpha }{2kk_1}}{\displaystyle \frac{\gamma _\alpha \widehat{k}\gamma \mu }{2kk_2}}\right],`$
$`\widehat{\mathrm{\Gamma }}_{\alpha \nu }`$ $`=`$ $`\left[\left({\displaystyle \frac{k_{1\alpha }}{kk_1}}{\displaystyle \frac{k_{1\alpha }}{kk_1}}\right)\gamma _\mu {\displaystyle \frac{\gamma _\alpha \widehat{k}\gamma \nu }{2kk_1}}{\displaystyle \frac{\gamma _\nu \widehat{k}\gamma \alpha }{2kk_2}}\right].`$ (28)
After contraction of the radiative leptonic tensor, we obtain

$$M_R^2=\frac{2e^6}{\stackrel{~}{Q}^4}L_{\mu \nu }^RW_{\mu \nu }=\frac{2e^6}{\stackrel{~}{Q}^4R_w}\underset{i=1}{\overset{5}{}}\theta _i_i.$$
(29)
$`R_w`$ was extracted explicitly in order to cancel $`v`$ coming from Jacobian (III). Arguments of $`_i`$ can be expressed as

$`\stackrel{~}{Q}^2`$ $`=`$ $`(qk)^2=Q^2+R_w\tau ,`$ (30)
$`\stackrel{~}{W}^2`$ $`=`$ $`(p+qk)^2=W^2R_w,`$
$`\stackrel{~}{t}`$ $`=`$ $`(qkp_h)^2=tR_w(\tau \mu ).`$ (31)
It is well known that the cross section of the radiative process is infrared–divergent, which requires careful consideration in order to cancel in the difference. The procedure will be discussed in the next section. Here we can extract the infrared convergent terms in $`\theta _i`$ in separate pieces:

$`\theta _i`$ $`=`$ $`{\displaystyle \frac{4}{R_w}}F_{IR}\theta _i^B+\theta _i^F.`$ (32)
The quantities $`\theta _i^B`$ are defined by expressions for $`\theta _i^0`$ with a reservation that the hadronic quantities $`V_1,V_2,S_t`$ and vector $`p_h`$ itself have to be calculated for the radiative kinematics (III). This term originates from the first terms of $`\mathrm{\Gamma }_{\mu \alpha }`$ and $`\overline{\mathrm{\Gamma }}_{\alpha \nu }`$ in definition of (III).

The explicit form of finite parts of these functions $`\theta _i^F`$ are given in Appendix.

Finally, the cross section of the radiative process ($`\sigma _R=d\sigma _R/dQ^2dW^2d\mathrm{\Omega }_h`$) is given as

$$\sigma _R=\frac{\alpha ^3}{2^9\pi ^4S^2W^4}𝑑\mathrm{\Omega }_k𝑑v\frac{\sqrt{\lambda _W}}{f}\underset{i}{}\frac{\theta _i_i}{\stackrel{~}{Q}^4},$$
(33)
where $`f=f_W/W`$.

As a cross–check, we consider the soft–photon limit $`\omega _{min}<R/2M<\omega _{max}`$ all energies and masses. In this case only the first term in (32) survives.

$$\sigma _R=\frac{2\alpha }{\pi }\left(\mathrm{log}\frac{Q^2}{m^2}1\right)\mathrm{log}\frac{w_{max}}{w_{min}}\sigma _0.$$
(34)
Integration over photonic angles is performed analytically,

$$𝑑\mathrm{\Omega }_kF_{IR}=2\left(l_m1\right),l_m=\mathrm{log}\frac{Q^2}{m^2}.$$
(35)
## IV Infrared divergence

As was already mentioned, this cross section contains an infrared divergence. Therefore, in order to compute RC, we have to use some regularization method first. We use the method of Bardin and Shumeiko for a covariant treatment of an infrared divergence problem. Basically, we follow the original papers devoted to this topic, namely, Ref. BSh  for $`e\mu `$ elastic scattering and Ref. Ak  where exclusive electroproduction was considered. One can find a good and detailed review in Ref. ABK .

Following the rules of dimensional regularization, we apply an identity transformation to the radiative cross section assuming that we deal with $`n`$-dimensional space:

$$\sigma _R=\sigma _R\sigma _{IR}+\sigma _{IR}=\sigma _F+\sigma _{IR}.$$
(36)
The purpose of the procedure is to separate the cross section into two pieces. The first term $`d\sigma _F`$ is complicated but infrared free. Infrared divergence is contained in the separate term $`d\sigma _{IR}`$, which has a quite simple structure and can be analytically calculated within $`n`$-dimensional space. In principle, there is some arbitrariness in choosing the form of the subtracted term. Actually only the asymptotical form for $`R_w`$ (or $`w`$) $`0`$ is fixed. Another limitation comes from the theorem about a possibility to switch the order of integration and the limit $`n4`$. It means that we have to provide uniform convergence of $`d\sigma _F`$ in the limit. Practically, it is ensured if the subtracted term has a structure $`F/w`$, and the difference is $`(F(w)F(0))/w`$.

We define the subtracted part of the radiative cross section as follows: Only the first term from Eq. (32) r.h.s. is kept. It gives the required $`1/R_w`$ behavior of the radiative cross section. Everywhere else, except for the $`\delta `$-function, we assume $`v=0`$. This allows to factorize the Born cross section when calculating the correction:

$$\sigma _{IR}=\sigma _0\frac{\alpha }{\pi ^2}𝑑v\frac{d\stackrel{}{k}}{\omega }F_{IR}\delta ((\mathrm{\Lambda }k)^2m_u^2).$$
(37)
In the cross section $`\sigma _F`$ we can now remove the regularization. The result is

$`\sigma _F`$ $`=`$ $`{\displaystyle \frac{\alpha ^3}{2^9\pi ^4S^2W^2}}{\displaystyle }d\mathrm{\Omega }_k{\displaystyle \frac{dv}{f}}\times `$
$`{\displaystyle \underset{i}{}}\left[{\displaystyle \frac{\sqrt{\lambda _W}}{\stackrel{~}{Q}^4}}\theta _i_i{\displaystyle \frac{4F_{IR}\sqrt{\lambda _W^0}}{Q^4}}\theta _i^0_i^0\right].`$
For calculation purposes, we split the integration region into two parts separated by the infinitesimal value of inelasticity $`\overline{v}`$:

$$\sigma _{IR}=\frac{\alpha }{\pi }\delta _R^{IR}\sigma _0=\frac{\alpha }{\pi }(\delta _S+\delta _H)\sigma _0$$
(39)
with

$`\delta _S`$ $`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle \underset{0}{\overset{\overline{v}}{}}}𝑑v{\displaystyle \frac{d^{n1}k}{(2\pi \mu )^{n4}k_0}F_{IR}\delta ((\mathrm{\Lambda }k)^2m_u^2)},`$
$`\delta _H`$ $`=`$ $`{\displaystyle \frac{1}{\pi }}{\displaystyle \underset{\overline{v}}{\overset{v_m}{}}}𝑑v{\displaystyle \frac{d^3k}{k_0}F_{IR}\delta ((\mathrm{\Lambda }k)^2m_u^2)}.`$ (40)
The term $`\delta _S^{IR}`$ corresponds to the soft photon contribution, while the term $`\delta _H^{IR}`$ is caused only by hard photons and therefore it does not contain the infrared singularity. Again, in the second contribution the regularization can be removed.

We keep integration over the 3–momentum of the radiated photon in a covariant form. It makes it possible to calculate the integral in any frame. We choose the frame where $`\stackrel{}{\mathrm{\Lambda }}=0`$ (so–called R–frame). Integration over inelasticity is external, therefore $`v`$ is fixed for the integral.

Calculation of the hard–photon contribution is straightforward,

$$\delta _H=(l_m1)\mathrm{log}\frac{v_m}{\overline{v}}.$$
(41)
For $`\delta _S`$, we follow Ref. BSh  generalizing these calculations to the exclusive electroproduction case.

Using the spherical $`n`$-dimensional frame, we have BSh ; Sirlin ; March

$`\delta _v={\displaystyle \frac{1}{\pi }}{\displaystyle \underset{0}{\overset{v_m}{}}}dv{\displaystyle \frac{2\pi ^{n/21}}{(2\pi \mu _0)^{n4}\mathrm{\Gamma }(\frac{n}{2}1)}}{\displaystyle }dk_0k_0^{n3}\times `$
$`\times {\displaystyle }\mathrm{sin}^{n3}\theta d\theta ({\displaystyle \frac{k_1}{2k_1k}}{\displaystyle \frac{k_2}{2k_2k}})^2\delta (v2k_0m_u).`$ (42)
$`F(\alpha ,cos(\theta ))=({\displaystyle \frac{k_1}{2k_1k}}{\displaystyle \frac{k_2}{2k_2k}})^2={\displaystyle \frac{1}{4k_0^2}}[{\displaystyle \frac{m^2}{E_{1R}^2(1\beta _1\mathrm{cos}\theta )^2}}`$
$`+{\displaystyle \frac{m^2}{E_{2R}^2(1\beta _2\mathrm{cos}\theta )^2}}{\displaystyle \underset{0}{\overset{1}{}}}{\displaystyle \frac{d\alpha Q^2}{E_\alpha ^2(1\beta _\alpha \mathrm{cos}\theta )^2}}]`$ (43)
Here Feynman parameters are introduced in order to join two denominators. Thus we define the new vector

$$k_\alpha =\alpha k_1+(1\alpha )k_2$$
(44)
and quantities $`\beta _{1,2,\alpha }`$ as a ratio of energy to momentum in the R–frame. In terms of invariants, we have

$`\beta _1=\sqrt{1+{\displaystyle \frac{4m^2m_u^2}{S_{0}^{}{}_{}{}^{2}}}},\beta _2=\sqrt{1+{\displaystyle \frac{4m^2m_u^2}{X_{0}^{}{}_{}{}^{2}}}},`$
$`\beta _\alpha =\sqrt{1+{\displaystyle \frac{4m_u^2(m^2+\alpha (1\alpha )Q^2)}{(\alpha S_0^{}+(1\alpha )X_0^{})^2}}}.`$ (45)
It should be noted that the same polar angles are used in Eq. (IV) for all three terms. That is why we can rotate the coordinate system for all the terms and choose the $`OZ`$ axis along $`\stackrel{}{k}_{1,2,\alpha }`$, respectively. After straightforward integration over the polar angle, we obtain Eq. (IV).

The next two steps involve integration over $`k_0`$ using the $`\delta `$-function and over $`v`$

$$\underset{0}{\overset{\overline{v}}{}}𝑑v𝑑k_0k_0^{n5}\delta (v2k_0m_u)=\frac{1}{n4}\left(\frac{\overline{v}}{2m_u}\right)^{n4}.$$
(46)
Since the pole $`(n4)^1`$ is extracted, we can use expansion in series of $`n4`$,

$`{\displaystyle \frac{1}{n4}}\left({\displaystyle \frac{\overline{v}}{4\pi \mu _0m_u}}\right)^{n4}{\displaystyle \frac{(\pi (1xi^2))^{n/22}}{\mathrm{\Gamma }(\frac{n}{2}1)}}`$ (47)
$`P_{IR}+\mathrm{log}{\displaystyle \frac{\overline{v}}{2\mu _0m_u}}+\mathrm{log}\sqrt{1\xi ^2}.`$
with $`\xi =cos\theta `$ and standard $`P_{IR}`$ defined as

$$P_{IR}=\frac{1}{n4}+\frac{1}{2}\gamma _E+\mathrm{ln}\frac{1}{2\sqrt{\pi }},$$
(48)
$`\gamma _E`$ being the Euler constant. After then we reduce the expression for $`\delta _S`$ to the form

$$\delta _S=𝑑\alpha 𝑑\xi \left[P_{IR}+\mathrm{log}\frac{\overline{v}}{2\mu _0m_u}+\mathrm{log}\sqrt{1\xi ^2}\right]F(\alpha ,\xi ),$$
(49)
which includes only standard integration that can be done using, for instance, tables from Appendix D of Ref.ABK . Integration is straightforward, but it contains an integral usually associated with a so–called $`S_\varphi `$ function in the Bardin and Shumeiko approach. Here we skip the discussion about properties of the function and give the result in the ultra–relativistic approximation which we use for the entire calculation:

$`S_\varphi `$ $`=`$ $`{\displaystyle \frac{Q^2}{2}}{\displaystyle \underset{0}{\overset{1}{}}}{\displaystyle \frac{d\alpha }{\beta _\alpha (m^2+\alpha (1\alpha )Q^2)}}\mathrm{log}{\displaystyle \frac{1\beta _\alpha }{1+\beta _\alpha }}=`$ (50)
$`{\displaystyle \frac{1}{2}}l_m^2l_m\mathrm{log}{\displaystyle \frac{S_0^{}X_0^{}}{m^2m_u^2}}{\displaystyle \frac{1}{2}}\mathrm{log}^2{\displaystyle \frac{S_0^{}}{X_0^{}}}`$
$`+\mathrm{Li}_2\left(1{\displaystyle \frac{Q^2m_u^2}{S_0^{}X_0^{}}}\right){\displaystyle \frac{\pi ^2}{3}}.`$
Finally, we have for $`\delta _S`$

$$\delta _S=2\left(P_{IR}+\mathrm{log}\frac{\overline{v}}{\mu M}\right)(l_m1)+\mathrm{log}\frac{S^{}X^{}}{m^2M^2}+S_\varphi .$$
(51)
The infrared–divergent terms $`P_{IR}`$, as well as the parameters $`\mu `$ and $`\overline{v}`$ are completely canceled in the sum $`\delta _S+\delta _H`$ with $`\delta _V`$ which is a contribution of the vertex function (Fig.1d):

$$\delta _V=2(P_{IR}+\mathrm{log}\frac{m}{\mu })(l_m1)\frac{1}{2}l_m^2+\frac{3}{2}l_m2+\frac{\pi ^2}{6}.$$
(52)
For this sum, we have

$$\frac{\alpha }{\pi }(\delta _S+\delta _H+\delta _V)=\delta _{inf}+\delta _{VR},$$
(53)
where

$`\delta _{VR}`$ $`=`$ $`{\displaystyle \frac{\alpha }{\pi }}\left({\displaystyle \frac{3}{2}}l_m2{\displaystyle \frac{1}{2}}\mathrm{ln}^2{\displaystyle \frac{X_0^{}}{S_0^{}}}+\mathrm{Li}_2\left[1{\displaystyle \frac{Q^2M^2}{S_0^{}X_0^{}}}\right]{\displaystyle \frac{\pi ^2}{6}}\right),`$
$`\delta _{inf}`$ $`=`$ $`{\displaystyle \frac{\alpha }{\pi }}(l_m1)\mathrm{ln}{\displaystyle \frac{v_m^2}{S_0^{}X_0^{}}}.`$ (54)
Higher–order corrections can be partially taken into account using a special procedure of exponentiation of multiple soft photon radiation. There is an uncertainty about which part of $`\delta _{VR}`$ has to be exponentiated. Within the approach Sh , we have to change $`1+\delta _{inf}`$ to $`\mathrm{exp}\delta _{inf}`$.

Collecting all the terms, we obtain for the case of meson electroproduction:

$$\sigma _{obs}=\sigma _0e^{\delta _{inf}}(1+\delta _{VR}+\delta _{vac})+\sigma _F.$$
(55)
Here the corrections $`\delta _{inf}`$ and $`\delta _{vac}`$ come from radiation of soft photons and effects of vacuum polarization, the correction $`\delta _{VR}`$ is an infrared–free sum of factorized parts of real and virtual photon radiation, and $`\sigma _F`$ is an infrared–free contribution from the bremsstrahlung process.

## V Leading Log Approximation

In this section we extract the leading log contribution from formulas obtained in the previous section. After then we show that the result coincides with what we obtain from a generalized leading log approximation.

The method of extraction of leading log (or peaking) contributions was first suggested in Ref. MoTsai  (see also papers Sh ; approx ). Formally, we have to calculate residues of the terms $`1/z_1`$ and $`1/z_2`$. The corresponding poles appear in $`\tau `$ for $`\tau =\tau _s=Q^2/S`$ and $`\tau =\tau _x=Q^2/X`$. As a result, we have (apart from the factorizable correction) two contributions of hard radiation,

$$\sigma _{LL}=\frac{\alpha }{2\pi }(l_m1)\left(\left(3+2\mathrm{log}\frac{v_m^2}{S^{}X^{}}\right)\sigma _0+\sigma _S+\sigma _X\right).$$
(56)
Technically, the contributions can be obtained as

$`{\displaystyle 𝑑\mathrm{\Omega }_k\theta _i}`$ $`=`$ $`8\pi l_m{\displaystyle \frac{W^2}{SQ^2}}{\displaystyle \frac{1+z_1^2}{z_1(1z_1)}}\theta _i^B,`$
$`{\displaystyle 𝑑\mathrm{\Omega }_k\theta _i}`$ $`=`$ $`8\pi l_m{\displaystyle \frac{W^2}{X+Q^2}}{\displaystyle \frac{1+z_2^2}{1z_2}}\theta _i^B.`$ (57)
For the two peaks the contributions are, respectively:

$`\sigma _s`$ $`=`$ $`{\displaystyle \underset{0}{\overset{v_m}{}}}{\displaystyle \frac{dv}{S^{}}}\left({\displaystyle \frac{1+z_1^2}{1z_1}}C_s\sigma _s{\displaystyle \frac{2\sigma _0}{1z_1}}\right),`$
$`\sigma _x`$ $`=`$ $`{\displaystyle \underset{0}{\overset{v_m}{}}}{\displaystyle \frac{z_2dv}{X^{}}}\left({\displaystyle \frac{1+z_2^2}{1z_2}}C_x\sigma _x{\displaystyle \frac{2\sigma _0}{1z_2}}\right),`$ (58)
where $`z_1=1v/S^{}`$, $`z_2=X^{}/(X^{}+v)`$ and

$`C_{s,x}={\displaystyle \frac{\sqrt{\lambda _W}}{\sqrt{\lambda _{Ws,x}^0}}}{\displaystyle \frac{W_{s,x}^2}{W^2}},`$
$`W_s^2=z_1SXz_1Q^2M^2,`$
$`W_s^2=S{\displaystyle \frac{X}{z_2}}{\displaystyle \frac{Q^2}{z_2}}M^2,`$ (59)
and $`\lambda _{Ws,x}^0`$ are calculated in accordance with Eq. (II) using $`W^2=W_{s,x}^2`$.

These leading log formulas were extracted from our exact expressions given in the previous section. The result can also be obtained using standard leading log techniques DeRujula  from our expression (III). The leading term comes from angular integration of denominators like $`kk_1`$ and $`kk_2`$. They can be extracted at the level of leptonic tensor

$$L_{\mu \nu }^R=L_{\mu \nu }\left[\frac{1+z_1^2}{1z_1}\frac{1}{kk_1}+\frac{1+z_2^2}{z_2(1z_2)}\frac{1}{kk_2}\right].$$
(60)
Only these scalar products in the denominators are subject to angular integration. For our exclusive process, integration yields:

$$\frac{d\stackrel{}{p}_u}{2ϵ_u}\frac{d\stackrel{}{k}}{2\omega }\delta (\mathrm{\Lambda }kp_u)\frac{1}{kk_{1,2}}=[\frac{\pi l_m}{S^{}},\frac{\pi l_m}{X^{}}].$$
(61)
Now, using (III), (60) and (61), we obtain the following leading–log cross section for the radiative process:

$`{\displaystyle \frac{d^6\sigma }{d\mathrm{\Gamma }^6}}`$ $`=`$ $`{\displaystyle \frac{\alpha l_m}{2\pi }}[{\displaystyle \frac{1}{S^{}}}{\displaystyle \frac{1+z_1^2}{1z_1}}{\displaystyle \frac{d^6\sigma _0(z_1k_1,k_2)}{d\stackrel{~}{\mathrm{\Gamma }}^6}}+`$
$`+{\displaystyle \frac{1}{X^{}}}{\displaystyle \frac{1+z_2^2}{1z_2}}{\displaystyle \frac{d^6\sigma _0(k_1,k_2/z_2)}{d\stackrel{~}{\mathrm{\Gamma }}^6}}],`$
$`d\mathrm{\Gamma }^6`$ $`=`$ $`{\displaystyle \frac{d^3\stackrel{}{k}_2}{2ϵ_2}}{\displaystyle \frac{d^3\stackrel{}{p}_h}{2ϵ_h}}.`$ (62)
While deriving the above results, (V,V), we take into account relations between $`z_{1,2}`$, $`v`$ and $`S^{}`$ ($`X^{}`$), which follow from the constraint $`(\mathrm{\Lambda }^2(1z_1)k_1)^2m_u^2=0`$ or $`(\mathrm{\Lambda }^2(11/z_2)k_2)^2m_u^2=0`$. They are

$$v=(1z_1)S^{},v=\frac{1z_2}{z_2}X^{}.$$
(63)
## VI Exclusive and inclusive radiative correction

An important consistency test is to show that inclusive RC can be obtained by integration over hadronic angles. It should be noted that this is not trivial because the hadronic angles of radiative and Born cross sections are defined even in different frames. For definition of the Born angles, we use the center-of-mass frame, while for angles of the radiative process we have to use another frame defined by vectors $`p`$ and $`qk`$. As a result, it leads to quite complicated kinematic relations between these angles.

We start with Eq.(56). After integration over hadronic angles, we have to obtain the leading log cross section ($`\sigma _{LL}^{inc}=d\sigma /dW^2dQ^2`$) for inclusive case. The inclusive formula can be found in Refs. approx ; POLRAD20 , for example. For double differential cross sections in $`Q^2`$ and $`W^2`$, it reads

$$\sigma _{LL}^{inc}=\frac{\alpha }{2\pi }(l_m1)((3+2\mathrm{log}\frac{v_m^2}{u_1u_2})\sigma _0^{inc}+\sigma _S^{inc}+\sigma _X^{inc},)$$
(64)
with

$`\sigma _S^{inc}`$ $`=`$ $`{\displaystyle \frac{\alpha }{2\pi }}l_m{\displaystyle \underset{z_1^m}{\overset{1}{}}}𝑑z_1z_1{\displaystyle \frac{1+z_1^2}{1z_1}}\sigma _{0S}^{inc}{\displaystyle \frac{2\sigma _0^{inc}}{1z_1}},`$ (65)
$`\sigma _X^{inc}`$ $`=`$ $`{\displaystyle \frac{\alpha }{2\pi }}l_m{\displaystyle \underset{z_2^m}{\overset{1}{}}}{\displaystyle \frac{dz_2}{z_2}}{\displaystyle \frac{1+z_2^2}{1z_2}}\sigma _{0X}^{inc}{\displaystyle \frac{2\sigma _0^{inc}}{1z_2}},`$ (66)
where $`z_1^m=X/u_2`$ and $`z_2^m=u_1/S`$.

Let us consider integration over $`\mathrm{\Omega }_h`$ of radiative cross sections ($`i.e.`$, the first terms on the left–hand–sides in (65)). The simplest way to relate angles in different frames is to express them in terms of kinematic invariants. We use the following relation,

$$d\mathrm{\Omega }_h=\frac{W^2}{\sqrt{\lambda _W^0}}\frac{dV_1^0dV_2^0}{\sqrt{D(V_1^0,V_2^0)}},$$
(67)
where $`D(V_1^0,V_2^0)=[ϵ(k_1,k_2,p,p_h)]^2`$ is Gramm determinant (see (14) and Ref.Bukling\_Kajanti ). For the cross sections $`\sigma _S^{inc}`$ and $`\sigma _X^{inc}`$ the corresponding expressions are

$$d\mathrm{\Omega }_h=\frac{W_s^2}{\sqrt{\lambda _{Ws}^0}}\frac{z_1dV_1dV_2}{\sqrt{D(z_1V_1,V_2)}}$$
(68)
and

$$d\mathrm{\Omega }_h=\frac{W_x^2}{\sqrt{\lambda _{Wx}^0}}\frac{dV_1dV_2}{z_2\sqrt{D(V_1,V_2/z_2)}}.$$
(69)
The expression for radiative kinematics is

$`16D(V_1,V_2)`$ $`=`$ $`(u_1V_2+u_2V_1(W^2+m_h^2m_u^2v)Q^2)^2`$ (70)
$`4(SXM^2Q^2)(V_1V_2m_h^2Q^2).`$
The Born case is reached in the limit $`v0`$. Two comments are in order before explicit integration. The integration area in variables $`V_1`$ and $`V_2`$ is defined by equation $`D(V_1,V_2)=0`$ that produces an ellipse. For radiative cross sections $`\sigma _S^{inc}`$ and $`\sigma _X^{inc}`$ this area is $`z_1`$ (or $`z_2`$)–dependent.

A simple integration procedure follows:

$`{\displaystyle 𝑑\mathrm{\Omega }_h\underset{0}{\overset{v_m}{}}𝑑v}{\displaystyle \underset{0}{\overset{v_m}{}}}𝑑v{\displaystyle 𝑑\mathrm{\Omega }_h}{\displaystyle \underset{0}{\overset{v_m}{}}}𝑑v{\displaystyle \frac{W^2}{\sqrt{\lambda _W}}}{\displaystyle \frac{dV_1dV_2}{\sqrt{D}}}`$
$`{\displaystyle 𝑑z_1𝑑V_1^s𝑑V_2\frac{W^2}{\sqrt{\lambda _W}}\frac{S^{}}{z_1\sqrt{D}}}`$
$`{\displaystyle 𝑑z_1𝑑\mathrm{\Omega }_h^s\frac{W^2}{\sqrt{\lambda _W}}\frac{S^{}}{z_1\sqrt{D}}\frac{\sqrt{\lambda _{Ws}}\sqrt{D_s}}{W_s^2}}.`$ (71)
The final formulas for the transformation

$$𝑑\mathrm{\Omega }_h\underset{0}{\overset{v_m}{}}𝑑v=\underset{z_1^m}{\overset{1}{}}𝑑z_1𝑑\mathrm{\Omega }_h^s\frac{S^{}}{C_s}$$
(72)
are obtained if we use the equality that can be checked directly:

$$\frac{D_s}{D}=z_1^2,\frac{D_x}{D}=\frac{1}{z_2^2}.$$
(73)
Using Eq. (72), one can see that the exclusive radiative cross section transforms to the inclusive one after integration.

Let us show how to obtain a factorized part of the inclusive correction from the exclusive result. The subtracted part of Eqs.(V) can be rewritten as

$`{\displaystyle \underset{0}{\overset{v_m}{}}}{\displaystyle \frac{dv}{S^{}(1z_1)}}`$ $`=`$ $`{\displaystyle \underset{z_{1^m}^1}{}}{\displaystyle \frac{dz_1}{1z_1}}{\displaystyle \underset{S^{}(v=v_m)}{\overset{S^{}(v=0)}{}}}{\displaystyle \frac{dS^{}}{S^{}}}`$ (74)
$`=`$ $`{\displaystyle \underset{z_{1^m}^1}{}}{\displaystyle \frac{dz_1}{1z_1}}\mathrm{log}\left({\displaystyle \frac{S_0^{}W}{u_1(Wm_h)}}\right).`$
After similar calculation for $`\sigma _X`$, we can see that inclusive RC is reproduced exactly.

## VII Numerical analysis

Based on the derived analytical formulas, a Fortran code EXCLURAD was developed. This code computes RC to the four–fold cross section ($`d^4\sigma /dWdQ^2d\mathrm{cos}\theta d\varphi _h`$) and to polarization beam asymmetry for processes (I). Both exact and leading log formulas obtained in the previous sections are included. Any value of the inelasticity cut can be optionally chosen.

Once the computational algorithm for RC is established, there are two possible ways of incorporating RC into analysis of experimental data on electroproduction:

a) Iteration procedure similar to the one implemented in the RC code POLRAD for inclusive reactions POLRAD20  or

b) Using realistic models for the structure functions (II) of coincidence electroproduction.

Although the former choice seems attractive due to its model independence, it requires, however, full and precise experimental mapping of the structure functions (II) in the entire range of the kinematic variables needed to compute the radiative process (33). Such a procedure is much more challenging and less efficient than for inclusive processes. Given that such mapping is not available at this time, the only choice left is (b). In this way, RC is applied to the model calculations, and then model parameters are fixed from available experimental data. Thus, RC appears as a necessary intermediate step in extracting model parameters from measurements.

We use the following models for calculation of structure functions (II):

MAID drec1 . In this model, baryon resonances are described using Breit–Wigner forms, while background contributions are described using standard Born terms, mixed pseudovector–pseudoscalar $`\pi `$NN coupling and t–channel vector meson exchange. The final amplitudes are constrained by unitarity and gauge invariance. We use two versions of MAID: the earlier one, denoted MAID’98 and the most recent, quoted here as MAID2000.

AO (Amplitude and Observables)burkert . The amplitudes are parameterized as follows: S–channel resonances are parameterized with relativistic Breit–Wigner forms with momentum–dependent widths. This part of the amplitude is complex. In addition, the s–channel and u–channel pion Born terms are included. These are real numbers. Additional real background amplitudes are used with an energy dependence that has correct threshold behavior.

Another model that can be included into the code EXCLURAD in a simple and straightforward way is the dynamical model of Sato and Lee sat01 . In this model, the off–shell non–resonant contributions to $`\gamma ^{}p\mathrm{\Delta }^+(1232)`$ were calculated directly by applying reaction theory within the Hamiltonian formulation underlying ’bare’ photo–coupling form factors. It should be noted that the Sato–Lee model does not include the contributions from higher resonances.

Since the computational algorithm for RC does not depend on a particular choice of hadronic structure functions, addition of any other model does not constitute a problem.

With the models at hand, we can now evaluate numerically the magnitude of RC as a function of various kinematics variables. Let us define the RC factor as follows,

$$\delta =\frac{\sigma _{obs}}{\sigma _0}.$$
(75)
For all the following plots, the electron beam energy is $`E_{beam}`$= 1.645 GeV and no cuts on inelasticity were used, except for Fig.9. We choose representative kinematics for current experiments with the CLAS detector at Jefferson Lab Hall B. The RC factors to the pion electroproduction cross sections are presented in Figs.3 and 4 as a function of $`W`$. In this region, the characteristic features of the cross section vs. $`W`$ are strong $`\mathrm{\Delta }`$(1232) and S<sub>11</sub>(1535) resonance peaks. One can see that RC is negative at the resonance peaks and positive in the dip regions between the peaks. Thus it tends to weaken the stength of the peaks and fill the dip regions above the given resonances. Also shown in Fig. 3 is the model dependence effect which appears to be noticeable at higher $`W`$ away from the resonance peaks. At higher values of $`W`$, the magnitude of RC is larger in the charged pion production case. This is mainly due to the wider range of inelasticity $`v`$ (24) associated with detection of a lighter hadron ($`i.e.`$, the pion).

Before discussing the angular dependence of RC, let us first comment on the definitions of hadronic angles in order to avoid possible ambiguities. The c.m.s. angles $`\theta `$ and $`\varphi `$ are between the direction of momentum lost by electrons ($`\stackrel{}{q}=\stackrel{}{k}_1\stackrel{}{k}_2`$) and momentum of the final pion $`\stackrel{}{p}_h`$, provided that the pion is detected. For the neutral pion production case, the convention is to also use the direction of pion momentum reconstructed from kinematics. In order to follow the convention for the neutral pion production, we define $`\theta `$ and $`\varphi `$ with respect to the direction of $`\stackrel{}{p}_h`$ opposite to the final proton momentum, while keeping in mind that it is different from the final pion momentum for the radiative processes Fig.1b,c.

The angular dependence of RC is shown in Figs.5 and 6, where it is plotted as a function of $`\mathrm{cos}\theta `$ and $`\varphi `$, respectively. The kinematics corresponds to the $`\mathrm{\Delta }`$(1232)-peak, where RC leads to suppression of the cross section. It can be seen from Fig.5 that, approaching the forward direction, at the given value of $`\varphi `$, the magnitude of RC as a function of $`\mathrm{cos}\theta `$ for the charged–pion channel smoothly increases from about 3% to 12%, while for the neutral pions is decreases from 10% to 5%. A common feature for both the channels is that RC is larger in magnitude for the parallel kinematics, when the detected hadron moves along the transferred momentum. The RC factor varies as a function of azimuthal angle $`\varphi `$, as well (Fig.6). Dependence of RC on the angle $`\varphi `$ has important implications for the super–Rosenbluth separation of electroproduction structure functions, and the $`\theta `$ dependence would affect the partial–wave analysis, resulting in corrections to electroexcitation parameters of baryon resonances.

The beam polarization asymmetry is plotted in Figs. 7 and 8 for the neutral and charged channels, respectively. One can see that for the asymmetry, RC changes from enhancement to suppression when passing across the $`\mathrm{\Delta }`$ and S<sub>11</sub> resonance regions. The RC factor is most substantial in the dip regions between resonances.

Figure 9 demonstrates the dependence on inelasticity cut $`v_{cut}`$ (24). It can be seen that for the smaller values of the cut, resulting in selection of softer bremsstrahlung photons, RC to both polarization–dependent and polarization–independent parts of the cross section is almost the same, resulting in a small correction to polarization asymmetry. As $`v_{cut}`$ increases, the correction to the asymmetry also increases due to the hard–photon emission coming into play.

Any experimental spectrum for electroproduction contains a radiative tail that, due to the finite energy resolution, cannot be experimentally separated from Born contribution. In practice, detected events are included into data analysis up to a certain (cut) limit and are interpreted as radiative events. The contribution of radiative events should be calculated theoretically using the same inelasticity cut ($`v_{cut}`$). The resulting contribution of the radiative events is subtracted from the integrated experimental spectrum in order to obtain the value of Born cross section that, naturally, should be cut–independent.

The Fortran code EXCLURAD can be downloaded from http://www.jlab.org/RC or obtained directly from the authors.

## VIII Discussion and conclusion

In this paper, we obtain explicit formulas for the lowest–order QED radiative correction to cross section and polarization beam asymmetry in the exclusive pion electroproduction. Analytic formulas are tested in several ways. Apart from traditional cross–checks like soft–photon and leading–log limits, it is found that integration with respect to the hadronic angles reproduces the inclusive radiative correction.

A Fortran code EXCLURAD is developed on the basis of the analytic formulas. Numerical analysis carried out for JLab kinematic conditions shows that:

* Radiative correction to the cross section of electroproduction is very sensitive to the cut on inelasticity. The harder cut eliminates contributions from the higher–energy part of the bremsstrahlung spectrum, thus leading to the smaller magnitude of RC to the polarization asymmetry.
* RC to cross sections can be as high as several tens of per cent.
* RC may vary depending on the chosen model for electroproduction structure functions. The proposed RC procedure may be viewed as a necessary intermediate step in interpretation of experimental results in terms of model parameters. Iteration procedure may be required for the regions where model dependence is substantial.
* RC have a non–trivial angular dependence in $`\mathrm{cos}\theta `$ and $`\varphi `$. This is particularly significant, as these are often used as input to extract structure functions and partial wave amplitudes.

## Acknowledgements

We would like to acknowledge useful discussions with our colleagues at Jefferson Lab. We thank Georgii Smirnov for reading the manuscript and useful comments. This work was supported by the US Department of Energy under contract DE-AC05-84ER40150.

## Appendix A

In this Appendix, we give formulas for $`\theta _i^F`$ ( Eq.(32)). In the unpolarized case,

$$\theta _i^F=\theta _{i2}+R_w\theta _{i3},$$
(76)
$`\theta _{12}`$ $`=`$ $`4F_{IR}\tau ,`$
$`\theta _{13}`$ $`=`$ $`4F2F_d,\tau ^2,`$
$`2\theta _{22}`$ $`=`$ $`4M^2F_{IR}\tau F_dS_p^2\tau +F_{1+}S_pS_x`$
$`+2F_2S_p+2F_{IR}S_x,`$
$`2\theta _{23}`$ $`=`$ $`4M^2F+2M^2F_d\tau ^2F_dS_x\tau F_{1+}S_p,`$
$`2\theta _{32}`$ $`=`$ $`4m_h^2F_{IR}\tau F_d\tau V_+^2+F_{1+}V_{}V_+`$
$`+2F_2\mu V_++2\mu F_{IR}V_{},`$
$`2\theta _{33}`$ $`=`$ $`4m_h^2F+2m_h^2F_d\tau ^2F_d\mu \tau V_{}F_{1+}\mu V_+,`$
$`2\theta _{42}`$ $`=`$ $`2F_dS_p\tau V_++F_{1+}S_pV_{}+F_{1+}S_xV_+`$
$`+2F_2\mu S_p+2F_2V_++2\mu F_{IR}S_x`$
$`4F_{IR}S_m\tau +2F_{IR}V_{},`$
$`2\theta _{43}`$ $`=`$ $`4FS_mF_d\mu S_x\tau +2F_dS_m\tau ^2F_d\tau V_{}`$
$`F_{1+}\mu S_pF_{1+}V_+.`$
Here $`V_\pm =V_1\pm V_2`$,

$`F_d`$ $`=`$ $`{\displaystyle \frac{1}{\kappa _1\kappa _2}},`$ (77)
$`F_{1+}`$ $`=`$ $`{\displaystyle \frac{1}{\kappa _1}}+{\displaystyle \frac{1}{\kappa _2}},`$ (78)
$`F_{2\pm }`$ $`=`$ $`\left({\displaystyle \frac{m^2}{\kappa _2^2}}\pm {\displaystyle \frac{m^2}{\kappa _1^2}}\right),`$ (79)
$`F_{IR}`$ $`=`$ $`F_{2+}Q^2F_d.`$ (80)
where $`\kappa _{1,2}`$ are defined in Eq. (27).

$`\theta _{52}`$ $`=`$ $`4(2F_{11}E_2+2F_{22}E_1+F_d(E_{12}\tau Q^2(E_1+E_2))`$
$`+{\displaystyle \frac{2}{S}}(E_1SE_{12}\tau +E_{12}E_2S)),`$
$`\theta _{53}`$ $`=`$ $`2(F_{1+}(E_2E_1)F_d(E_1+E_2)\tau `$ (81)
$`+{\displaystyle \frac{4}{S}}F_{11}E_2(\tau 1)),`$
$`E_{12}`$ $`=`$ $`ϵ(k_1,k_2,p_1,p_h)={\displaystyle \frac{1}{4}}[(V_1XQ^2S_m+SV_2)^2+`$
$`+4(SXM^2Q^2)(Q^2m_h^2V_1V_2)]^{1/2},`$
$`E_1`$ $`=`$ $`ϵ(k,k_1,p_1,p_h)={\displaystyle \frac{1}{4}}[(R_wV_1(\tau 1)R_w\mu S+S_m\kappa _1)^2+`$
$`+4(\kappa _1M^2+R_w(\tau 1)S)(\kappa _1m_h^2+V_1R_w\mu )]^{1/2},`$
$`E_2`$ $`=`$ $`ϵ(k,k_2,p_1,p_h)={\displaystyle \frac{1}{4}}[(R_wV_2(\tau 1)R_w\mu X+S_m\kappa _2)^2+`$ (82)
$`+4(\kappa _2M^2+R_w(\tau 1)X)(\kappa _2m_h^2+V_2R_w\mu )]^{1/2}.`$


text:
# Projected BCS Wave Functions for Low Dimensional Frustrated Spin Systems

11institutetext: Kavli Institute for Theoretical Physics and Department of Physics, University of California, Santa Barbara, CA 9310611institutetext: Istituto Nazionale per la Fisica della Materia and Université de Lausanne, CH-1015 Lausanne, Switzerland11institutetext: Istituto Nazionale per la Fisica della Materia and Università dell’Insubria, I-22100 Como, Italy11institutetext: Istituto Nazionale per la Fisica della Materia, and SISSA, I-34014 Trieste, Italy

Non-magnetic ground states are a fascinating possibility allowed by the physics of quantum antiferromagnets. These states – which lack the classical periodical long-range order – can be stabilized whenever reduced dimensionality, a small spin value, and/or the presence of competing interactions lead to strong enough zero-point quantum fluctuations.

One dimensional or quasi-one dimensional spin-half Heisenberg antiferromagnets often have non-magnetic ground states. Indeed for the 1D nearest-neighbors Heisenberg model

$$\widehat{}=J{}_{1}{}^{}\underset{n.n.}{}\widehat{𝐒}_i\widehat{𝐒}_j$$
(1)
the exact solution due to Bethe  predicts the absence of true long-range antiferromagnetic order even if the ground state is very close to have a broken symmetry, with a gapless excitation spectrum and a power-law decay of spin-spin correlations. This is also the case for any array consisting of an odd number of chains (odd-leg ladder systems). The ground state of two chains or in general of any even-leg ladder system is non-magnetic too. However – in contrast to the previous cases – here the correlation length is finite and the spectrum is gapped. Such a gap is known to decrease exponentially with the number of legs  leading to a gapless spectrum in the two dimensional limit, where the ground state of the Heisenberg model has genuine long-range antiferromagnetic order.

Competing interactions may allow in principle the stabilization of a non-magnetic ground state even in truly two dimensional systems. One of the simplest examples of these frustrated systems, which has been also recently realized experimentally , is the so-called $`J_1J_2`$ model

$$\widehat{}=J{}_{1}{}^{}\underset{n.n.}{}\widehat{𝐒}_i\widehat{𝐒}_j+J{}_{2}{}^{}\underset{n.n.n.}{}\widehat{𝐒}_i\widehat{𝐒}_j,$$
(2)
where the antiferromagnetic alignment between neighboring spins (due to $`J_1>0`$) is hindered by a next-nearest-neighbors antiferromagnetic coupling ($`J_2>0`$).

Classically, the minimum energy configuration of the 2D $`J_1J_2`$ model has the conventional Néel order with magnetic wave vector $`𝐐=(\pi ,\pi )`$ for $`J_2/J_1<0.5`$. Instead for $`J_2/J_1>0.5`$ the minimum energy configuration is the so-called collinear state with the spins ferromagnetically aligned in one direction and antiferromagnetically in the other, corresponding to magnetic wave vectors $`𝐐=(\pi ,0)`$ or $`𝐐=(0,\pi )`$. Exactly at $`J_2/J_1=0.5`$ any classical state having zero total spin on each elementary square plaquette is a minimum of the total energy. These states include both the Néel and the collinear states but also many others with no long-range order so that the occurrence of a non-magnetic ground state in the quantum case, for a small spin value, is likely around this value of the $`J_2/J_1`$ ratio. Indeed, at present there is a general consensus on the fact that the combined effect of frustration and zero-point motion leads to the disappearance of the long-range antiferromagnetic order marked by the opening of a finite spin gap for $`0.4<J_2/J_1<0.55`$. The nature of this non-magnetic ground state is one of the most interesting puzzles of the physics of frustrated spin systems. In particular an open question is whether the ground state of the $`J_1J_2`$ Heisenberg model is a homogeneous spin liquid, i.e., a state with all the symmetries of the Hamiltonian, as it was originally suggested by Figueirido et al. . The other possibility is a ground state which is still $`𝒮𝒰(2)`$ invariant, but nonetheless breaks some crystal symmetries, dimerizing in some special pattern (see below).

A simple picture of a non-magnetic ground state can be given in terms of the so-called Resonating Valence Bond (RVB) states. These are linear superpositions of valence bond states in which each spin forms a singlet bond with another spin on the opposite sublattice (say $`A`$ and $`B`$)

$$|\psi _{RVB}=\underset{i_\alpha A,j_\beta B}{}h(r_1)\mathrm{}h(r_{\frac{N}{2}})(i_1j_1)\mathrm{}(i_{\frac{N}{2}},j_{\frac{N}{2}}),$$
(3)
where $`N`$ is the number of sites of the lattice, $`r_m`$ is the distance between the spins forming the $`m^{th}`$ singlet bond $`(i_mj_m)`$, and $`h(r_m)`$ is a bond weight factor. These states form in general a (overcomplete) basis of the $`S=0`$ subspace so that any singlet wave function can be represented in terms of them. However, they represent a non-magnetic state whenever the short-ranged bonds dominate the superposition (3). More precisely, it has been numerically shown by Liang, Doucot and Anderson that the RVB state (3) has no long-range antiferromagnetic order for bonds that decay as rapidly as $`h(r)r^p`$, with $`p5`$. Such bonds can be either homogeneously spatially distributed on the lattice, with short-range correlations among each other (spin liquid) \[fig. 1 (a)\], or they can break some symmetries of the Hamiltonian, with the dimers frozen in some special patterns \[fig. 1 (b)\] as originally predicted for the $`J_1J_2`$ model in the regime of strong frustration.

In a seminal paper, Anderson proposed that a physically transparent description of a RVB state can be obtained in fermionic representation by starting from a BCS-type pairing wave function, of the form

$$|\psi _{\mathrm{BCS}}=\mathrm{exp}\left(\underset{i,j}{}f_{i,j}\widehat{c}_{i,}^{}\widehat{c}_{j,}^{}\right)|0.$$
(4)
This wave function is the ground state of the well-known BCS Hamiltonian

$$H_{BCS}=\underset{k,\sigma }{}ϵ_kc_{k,\sigma }^{}c_{k,\sigma }+\underset{k}{}(\mathrm{\Delta }_kc_{k,}^{}c_{k,}^{}+h.c.)$$
(5)
where $`\sigma =,`$, $`ϵ_k=2[\mathrm{cos}k_x+\mathrm{cos}k_y]`$ is the free-electron dispersion and $`\mathrm{\Delta }_k=\mathrm{\Delta }_k`$ is the (real) gap function, provided the Fourier transform $`f_k`$ of the pairing function, $`f_{i,j}`$, satisfies: $`f_k=\mathrm{\Delta }_k/(ϵ_k+\sqrt{ϵ_k^2+\mathrm{\Delta }_k^2})`$. The non-trivial character of this wave function emerges when we restrict to the subspace of fixed number of electrons (equal to the number of sites) and enforce Gutzwiller projection onto the subspace with no double occupancies: singlet pairs do not overlap in real space and this wave function can be described by a superposition of valence bond states of the form (3).

This projected-BCS (p-BCS) wave function turns out to be an almost exact representation of several low-dimensional spin system with non-magnetic ground states. For instance it provides an excellent variational ansatz of the ground state of the Heisenberg chain and of the two-leg ladder giving a very accurate estimate of the ground-state energy (fig. 2) and reproducing almost exactly the antiferromagnetic correlations. In the first case the spin structure factor $`S(q)=\widehat{𝐒}_q\widehat{𝐒}_q`$ shows a cusp at $`q=\pi `$ while for two-leg ladders it has a broad maximum at $`q=(\pi ,\pi )`$. These features are remarkably well reproduced by the (p-BCS) variational wave function (fig. 2), which generates robust antiferromagnetic correlations at short distances with a very simple parameterization of the gap function: $`\mathrm{\Delta }_k=\mathrm{\Delta }_1\mathrm{cos}k+\mathrm{\Delta }_2\mathrm{cos}3k`$ for the chain and $`\mathrm{\Delta }_k=\mathrm{\Delta }_x\mathrm{cos}k_x+\mathrm{\Delta }_y\mathrm{cos}k_y`$ for the ladder.

In two dimensions this wave function has been already studied for the pure Heisenberg model by several authors for $`\mathrm{\Delta }_k(\mathrm{cos}k_x\mathrm{cos}k_y)`$. In this case it provides a reasonable prediction for the ground-state energy but it fails in reproducing correctly the long-range antiferromagnetic order of the ground state. Here we show that this type of RVB state represents an extremely accurate variational ansatz for the $`J_1J_2`$ model in the non-magnetic phase when the gap function $`\mathrm{\Delta }_k`$ is carefully parameterized. In particular, a definite symmetry is guaranteed to the p-BCS state provided $`\mathrm{\Delta }_k`$ transforms according to a one dimensional representation of the spatial symmetry group. A careful analysis  shows that the odd component of the gap function $`\mathrm{\Delta }_k=\mathrm{\Delta }_{k+(\pi ,\pi )}`$ may have spatial symmetries different from those of the even component $`\mathrm{\Delta }_k=\mathrm{\Delta }_{k+(\pi ,\pi )}`$. Indeed, the best variational energy is obtained when the former has $`d_{x^2y^2}`$ symmetry, whereas the latter either vanishes or it has $`d_{xy}`$ symmetry. In order to determine the best variational wave function of this form we have used a recently developed quantum Monte Carlo technique  that allows to optimize a large number of variational parameters with modest computational effort.

The remarkable accuracy of the p-BCS wave function in describing the ground state of the 2D $`J_1J_2`$ model in the regime of strong frustration can be shown by calculating the variational energy and the overlap with the exact ground state, $`|\psi _0`$, for the largest square cluster $`N=6\times 6`$ where the solution can be numerically determined by exact diagonalization. As shown in fig. 3 the accuracy of the p-BCS wave function rapidly increase by increasing the frustration ratio $`J_2/J_1`$ whereas conventional Néel ordered spin-wave wave functions  quickly become less and less accurate. Entering the regime of strong frustration $`J_2/J_10.45\pm 0.05`$, where a gapped non-magnetic ground state is expected, the p-BCS wave function becomes impressively accurate with a relative accuracy on the ground-state energy of order $`4\times 10^3`$ and an overlap to the exact ground state of $`99\%`$, both improved by more than an order of magnitude with respect to the $`J_2=0`$ case. This fact implies that the ground state in the strongly frustrated regime is almost exactly reproduced by a RVB wave function.

Interestingly, the transition to the regime of strong frustration is marked by the stabilization at the variational level of a non-zero $`d_{xy}`$ component of the gap function. This allows to reproduce correctly the phases of the actual ground-state configurations as illustrated in fig. 4. A measure of the accuracy of a variational wave function $`|\psi _V`$ in reproducing the phases of the ground state can be given in terms of the average sign $`S=_x|x|\psi _\mathrm{V}|^2\mathrm{Sgn}\left[x|\psi _\mathrm{V}x|\psi _0\right].`$ In the unfrustrated case it is well known that such phases are determined by the so-called Marshall sign rule : on each real space configuration $`|x`$, the sign of the ground-state wave function is determined only by the number of spin down in one of the two sublattices. This feature, rigorously valid for $`J_2=0`$, turns out to be a very robust property for weak frustration ($`J_2/J_1<0.3`$). However, it is clearly violated when frustration plays an important role. It can be shown  that the Marshall sign (i.e., $`S=1`$ for $`J_2/J_1=0`$) can be obtained using the p-BCS wave function, with only the $`d_{x^2y^2}`$ component, so that this wave function provides an almost exact representation of the ground-state phases for weak frustration. However, for $`J_2/J_1>0.4`$, the phases of the wave function are considerably affected by the strong frustration and only when a sizable $`d_{xy}`$ component is stabilized at the variational level, this property can be correctly reproduced.

An even more remarkable effect associated with the $`d_{xy}`$ component of the gap function is the change induced on antiferromagnetic correlations. As it is shown, in fig. 5 the finite-size magnetic structure factor of the $`d_{x^2y^2}`$ p-BCS state is sharply peaked around the antiferromagnetic wave vector $`𝐐=(\pi ,\pi )`$ giving rise to a logarithmic divergence in the thermodynamic limit. Such a divergence is instead washed out in presence of the $`d_{xy}`$ component of the gap function, leading to a state with weaker short-range antiferromagnetic correlations which is of course more suitable to describe the spin-gapped strongly-frustrated phase.

Of course, the accuracy in the energy does not necessarily guarantee a corresponding accuracy in correlation functions. However, as shown in fig. 6, the comparison of the magnetic structure factor with the exact result gives a clear indication that correlation functions obtained by the variational approach are essentially exact. Furthermore using the stochastically implemented Lanczos technique and the variance extrapolation method  we have verified that the accuracy of the p-BCS wave function is preserved by increasing the lattice size.

In order to investigate the existence of long-range dimer-like correlations, as in the columnar or the plaquette valence bond state, we have calculated the dimer-dimer correlation functions, $`\mathrm{\Delta }_{i,j}^{k,l}=\widehat{S}_i^z\widehat{S}_j^z\widehat{S}_k^z\widehat{S}_l^z\widehat{S}_i^z\widehat{S}_j^z\widehat{S}_k^z\widehat{S}_l^z`$. In presence of some broken spatial symmetry, the latter should converge to a finite value for large distance. This is clearly ruled out by our results, shown in fig. 7, with a very robust confirmation of the liquid character of the ground state for $`J_2/J_10.5`$, which is correctly described by our variational approach.

A totally symmetric spin-liquid solution proposed for this model in ref.  was actually rather unexpected after the work of Read and Sachdev, providing arguments in favor of spontaneous dimerization. This conclusion was supported by series expansion  and quantum Monte Carlo studies included the one done by two of us. It is clear however that it is very hard to reproduce a fully symmetric spin liquid ground state, with any technique, numerical or analytical, based on reference states explicitly breaking some lattice symmetry.

In conclusion, the spin-liquid RVB ground state, originally proposed to explain high-Temperature superconductivity, is indeed a very robust property of strongly frustrated low-dimensional spin systems. Due to the success in reproducing the non-magnetic ground states of other low-dimensional spin systems like the 1D chain and the two-leg ladder,  we expect that the p-BCS RVB wave function represents the generic variational state for a spin-half spin liquid, once the pairing function $`f_{i,j}`$ is exhaustively parameterized according to the symmetries of the Hamiltonian. Work is in progress on this line of research.

###### Acknowledgements.

We thank C. Lhuillier, F. Mila, and D.J. Scalapino for stimulating discussions. This work has been partially supported by MURST (COFIN01). L.C. was supported by NSF grant DMR-9817242.