format: text

text:
# I Introduction

## I Introduction

The deeply virtual Compton scattering (DVCS) process, in which a highly virtual photon $`\gamma ^{}(q_1)`$ produces a real photon $`\gamma (q_2=q_1+r)`$ at small invariant momentum transfer $`t=r^2`$, is receiving a lot of attention as a potential source of new information about the nucleon structure in terms of so-called “skewed” parton distributions (SPD’s) . The first experimental observation of DVCS has been already reported by the ZEUS collaboration . Originally, only the purely twist–2 contribution to the tensor amplitude $`T_{\mu \nu }`$ amplitude was calculated , which does not include subleading terms proportional to the transverse part $`r_{}\mathrm{\Delta }`$ of the momentum transfer. As a result, the twist-2 amplitude formally violates the electromagnetic (EM) gauge invariance in terms linear in $`\mathrm{\Delta }`$. Guichon and Vanderhaeghen  proposed to add an ad hoc $`O(\mathrm{\Delta })`$ term producing an “improved” expression which is EM gauge invariant up to $`O(\mathrm{\Delta }^2)`$ terms. Their prescription was recently supported by several groups  which derived this term in a regular way as a kinematical twist-3 contribution. Furthermore, Anikin et al. , using the momentum–space collinear expansion, have obtained expressions which include (in case of the pion target) all the relevant twist–3 operators. The DVCS amplitude to accuracy $`1/\sqrt{q_1^2}`$ was also calculated by Penttinen et al.  in a parton model approach. Within the light-cone expansion framework, Belitsky and Müller  analyzed both quark and quark-gluon contributions and demonstrated that, to get a gauge invariant result up to terms of order $`t/q_1^2`$, it is sufficient to retain only the part of the twist–3 SPD’s which is obtained by Wandzura–Wilczek (WW) type formulas from the twist–2 distributions. Recently, Kivel et al.  established that the WW-approximation expression for the transverse polarization of the virtual photon diverges. The mathematical aspects of twist decomposition were discussed by Blümlein and Robaschik , Geyer and Lazar .

Our goal is to analyze the DVCS amplitude within the QCD string operator approach of Balitsky and Braun , which proved to be a powerful tool to investigate the higher-twist effects. Here, we consider only the kinematical twist-3 terms. In this sense, our results are equivalent to the WW approximation. In addition to offering an alternative derivation of this approximation, we incorporate the formalism of double distributions  which provides a simple way of deriving relations between new SPD’s. For instance, the fact that the WW approximation for the amplitude $`T_{\mu \nu }`$ gives finite results for longitudinally polarized photon, but diverges in the case of transverse polarization, can be easily understood on the basis of the formulas relating these SPD’s and the basic twist-2 DD’s. We also note that the divergent part of $`T_{\mu \nu }`$ has a vanishing projection onto the polarization of the final real photon.

## II DVCS amplitude

Generalities. The virtual Compton scattering amplitude is derived from the correlation function

$`T_{\mu \nu }`$ $`=`$ $`i{\displaystyle d^4xd^4ye^{i(q_1x)+i(q_2y)}p_2|\mathrm{T}\left\{J_\mu (x)J_\nu (y)\right\}|p_1},`$ (1)
where $`J^\mu (x)`$ is the electromagnetic current operator. Due to current conservation $`^\mu J_\mu (x)=0`$, this function is transverse with respect to the incoming and outgoing photon momenta: $`q_1^\mu T_{\mu \nu }=0,q_2^\nu T_{\mu \nu }=0.`$ It is convenient to switch to symmetric variables $`q=(q_1+q_2)/2`$, $`r=q_2q_1`$ and $`p=(p_1+p_2)/2`$. Then, the transversality conditions convert into two relations

$`q^\mu T_{\{\mu \nu \}}={\displaystyle \frac{r^\mu }{2}}T_{[\mu \nu ]},q^\mu T_{[\mu \nu ]}={\displaystyle \frac{r^\mu }{2}}T_{\{\mu \nu \}}`$ (2)
connecting the symmetric $`T_{\{\mu \nu \}}(T_{\mu \nu }+T_{\nu \mu })/2`$ and antisymmetric $`T_{[\mu \nu ]}(T_{\mu \nu }T_{\nu \mu })/2`$ parts of $`T^{\mu \nu }`$. In the $`r=0`$ forward limit, the two relations decouple to give the DIS transversality conditions $`q^\mu T_{\{\mu \nu \}}=0`$, $`q^\mu T_{[\mu \nu ]}=0`$.

In DVCS, the initial photon is in the Bjorken kinematics {$`q_1^2\mathrm{},(p_1q_1)\mathrm{},x_Bq_1^2/[2(p_1q_1)]`$ fixed} and the final one is real $`q_2^2=0`$. Since $`q_2^2=q_1^2+2(q_1r)+t`$, the momentum transfer $`r`$ in this process should have a large component in the direction of $`p`$, with $`(rq_1)`$ close to $`q_1^2/2`$ for small $`t`$. The size of this component is characterized by the skewedness parameter $`\eta (rq)/2(pq)`$. For DVCS, $`\eta `$ coincides, up to $`O(t)`$ terms, with the generalized Bjorken variable $`\xi q^2/2(pq).`$ Hence, the momentum transfer may be split into the component parallel to $`p`$ and the remainder $`\mathrm{\Delta }`$

$`r=2\xi p+\mathrm{\Delta },`$ (3)
which in the $`t=0`$, $`p^2=0`$ limit is transverse both to $`p`$ and $`q`$: $`(\mathrm{\Delta }q)=t/4`$, $`(\mathrm{\Delta }p)=2\xi p^2`$.

Coordinate representation. An efficient way to study the behavior of Compton amplitudes in the Bjorken limit is to use the light–cone expansion for the product $`\mathrm{\Pi }_{\mu \nu }(x,y)i\mathrm{T}J_\mu (x)J_\nu (y)`$ of two vector currents in the coordinate representation. Following Balitsky and Braun , we start from the formal light–cone expansion in terms of QCD string operators (with gauge links along the straight line between the fields which, for brevity, we do not write explicitly). The leading light–cone singularity is contained in the “handbag” contribution

$`\mathrm{\Pi }_{\mu \nu }(z|X)={\displaystyle \frac{4iz_\rho }{\pi ^2z^4}}\{s_{\mu \rho \nu \sigma }𝒪_\sigma (z|X)ϵ_{\mu \rho \nu \sigma }𝒪_{5\sigma }(z|X)\},`$ (4)
where $`s_{\mu \rho \nu \sigma }=g_{\mu \rho }g_{\nu \sigma }g_{\mu \nu }g_{\rho \sigma }+g_{\mu \sigma }g_{\nu \rho }`$, $`X=(x+y)/2`$, $`z=yx`$, and

$`𝒪_\sigma (z|X)`$ $`=`$ $`{\displaystyle \frac{1}{2i}}\left[\overline{\psi }(Xz/2)\gamma _\sigma \psi (X+z/2)(zz)\right],`$ (5)
$`𝒪_{5\sigma }(z|X)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\overline{\psi }(Xz/2)\gamma _\sigma \gamma _5\psi (X+z/2)+(zz)\right].`$ (6)
Twist–2 part. The string operators in Eq. (4) do not have a definite twist. The twist–2 part is defined by formally Taylor–expanding the string operators in Eq. (4) in the relative coordinate $`z`$ and retaining only the totally symmetric traceless parts of the coefficients in the expansion:

$`\left[\overline{\psi }(Xz/2)\gamma _\sigma \psi (X+z/2)\right]^{\mathrm{twist}2}{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{1}{n!}}z_{\mu _1}\mathrm{}z_{\mu _n}\overline{\psi }(X)\left[\gamma _{\{\sigma }\underset{\mu _1}{\overset{}{D}}\mathrm{}\underset{\mu _n\}}{\overset{}{D}}\text{traces}\right]\psi (X),`$ (7)
and similarly for the operator with $`\gamma _\sigma \gamma _5`$. As shown in Ref., “symmetrization” and “subtraction of traces” can be carried out directly at the level of non-local operators. The part of the string operator corresponding to totally symmetric local tensor operators is projected out by

$$\left[\overline{\psi }(Xz/2)\gamma _\sigma \psi (X+z/2)\right]^{\mathrm{sym}}=\frac{}{z_\sigma }_0^1𝑑t\overline{\psi }(Xtz/2)\widehat{z}\psi (X+tz/2)$$
(8)
(we use the notation $`\widehat{z}z^\sigma \gamma _\sigma `$). The subtraction of traces in the local operators implies that the twist–2 string operator contracted with $`z_\sigma `$ should satisfy the d’Alembert equation with respect to $`z`$:

$`\mathrm{}_z\left[\overline{\psi }(Xtz/2)\widehat{z}\psi (X+tz/2)\right]^{\mathrm{twist}2}=0.`$ (9)
Transversality and twist–3 operators.  In the coordinates $`X`$ and $`z`$, the transversality conditions (2) are

$`{\displaystyle \frac{}{z_\mu }}\mathrm{\Pi }_{\{\mu \nu \}}(z|X)={\displaystyle \frac{1}{2}}{\displaystyle \frac{}{X_\mu }}\mathrm{\Pi }_{[\mu \nu ]}(z|X),{\displaystyle \frac{}{z_\mu }}\mathrm{\Pi }_{[\mu \nu ]}(z|X)={\displaystyle \frac{1}{2}}{\displaystyle \frac{}{X_\mu }}\mathrm{\Pi }_{\{\mu \nu \}}(z|X).`$ (10)
Consider the part of the current product given by Eq. (4) with the string operators replaced by their twist–2 parts. From Eq. (9) and $`(/z_\rho )[z_\rho /(2\pi ^2z^4)]=i\delta ^{(4)}(z)`$ it follows that $`(/z_\mu )\mathrm{\Pi }_{\{\mu \nu \}}^{\mathrm{twist}2}=0`$, $`(/z_\mu )\mathrm{\Pi }_{[\mu \nu ]}^{\mathrm{twist}2}=0.`$ Since forward matrix elements are zero for all total derivative operators, this guarantees the transversality of the twist–2 contribution in the case of deep inelastic scattering. In the non-forward case, we have $`(/X_\mu )\mathrm{\Pi }_{\{\mu \nu \}}^{\mathrm{twist}2}0,(/X_\mu )\mathrm{\Pi }_{[\mu \nu ]}^{\mathrm{twist}2}0,`$ and (10) is violated. The non-transverse terms in the twist–2 contribution can only be compensated by contributions from operators of higher twist. In fact, the necessary operators are contained in the part of the string operator which was dropped in taking the twist–2 part. Incorporating QCD equations of motion, it is possible to show  that the twist$`>2`$ part involves the total derivatives of string operators

$`\overline{\psi }(z/2)\gamma _\alpha \psi (z/2)\left[\overline{\psi }(z/2)\gamma _\alpha \psi (z/2)\right]^{\mathrm{sym}}={\displaystyle \frac{i}{2}}ϵ_{\alpha \xi \rho \kappa }z_\xi {\displaystyle \frac{}{X_\rho }}{\displaystyle _0^1}𝑑tt\overline{\psi }(tz/2)\gamma _\kappa \gamma _5\psi (tz/2)+\mathrm{}.`$ (11)
The ellipses stand for quark–gluon operators (we do not write them explicitly since they are not needed to restore transversality of the twist–2 contribution). The relation for the operator with Dirac matrix $`\gamma _\alpha \gamma _5`$ is obtained by changing $`\gamma _\alpha \gamma _\alpha \gamma _5`$, $`\gamma _\kappa \gamma _5\gamma _\kappa `$. The operators appearing under the total derivative on the R.H.S. of Eq. (11) and its $`\gamma _\alpha \gamma _5`$ analog are still the full string operators with no definite twist. Hence, one can decompose them into a symmetric (i.e., twist–2) part and total derivatives, and so on; thus expressing the original string operator as the sum of its symmetric part and an infinite series of arbitrary order total derivatives of symmetric operators. This series can be summed up in a closed form (the details of our calculation are presented elsewhere ; similar expressions were derived independently in ). Up to operators whose matrix elements give $`O(t)`$ contributions to the Compton amplitude, the result is

$`\overline{\psi }(z/2)\gamma _\sigma \psi (z/2)`$ $`=`$ $`{\displaystyle _0^1}𝑑v\left\{\mathrm{cos}\left[{\displaystyle \frac{i\overline{v}}{2}}\left(z{\displaystyle \frac{}{X}}\right)\right]{\displaystyle \frac{}{z_\sigma }}+{\displaystyle \frac{iv}{2}}\mathrm{sin}\left[{\displaystyle \frac{i\overline{v}}{2}}\left(z{\displaystyle \frac{}{X}}\right)\right]{\displaystyle \frac{}{X_\sigma }}\right\}\overline{\psi }(vz/2)\widehat{z}\psi (vz/2)`$ (13)
$`+{\displaystyle \frac{i}{2}}ϵ_{\sigma \alpha \beta \gamma }z_\alpha {\displaystyle \frac{}{X_\beta }}{\displaystyle \frac{}{z_\gamma }}{\displaystyle _0^1}𝑑v{\displaystyle _v^1}𝑑u\mathrm{cos}\left[{\displaystyle \frac{i\overline{u}}{2}}\left(z{\displaystyle \frac{}{X}}\right)\right]\overline{\psi }(vz/2)\widehat{z}\gamma _5\psi (vz/2)+\mathrm{}.`$
An analogous formula applies to the operators with $`\gamma _\sigma \gamma _\sigma \gamma _5`$; one should just replace $`\widehat{z}\widehat{z}\gamma _5,\widehat{z}\gamma _5\widehat{z}`$.

## III Parametrization of nonforward matrix elements

Double distributions. To get the amplitude for deeply virtual Compton scattering off a hadronic target we need parametrizations of the hadronic matrix elements of the uncontracted twist–2 string operators $`𝒪_\sigma ,𝒪_{5\sigma }`$ appearing in Eq. (4). We will derive them from Eq. (13). For simplicity, we consider here one quark flavor and the pion target, which has zero spin and practically vanishing mass. In this case, the matrix element of the contracted axial operator $`z^\sigma 𝒪_{5\sigma }(z|\mathrm{\hspace{0.17em}0})`$ (parametrized in the forward limit by the polarized parton density) is identically zero. Thus we need only the parametrization for the matrix element $`pr/2|𝒪(z|\mathrm{\hspace{0.17em}0})|p+r/2`$ of the contracted vector operator $`𝒪(z|\mathrm{\hspace{0.17em}0})z^\sigma 𝒪_\sigma (z|\mathrm{\hspace{0.17em}0})`$. With respect to $`z`$, it can be regarded as a function of three invariants $`(pz),(rz)`$ and $`z^2`$. For dimensional reasons, the dependence on $`z^2`$ is through the combinations $`tz^2`$ and $`p^2z^2`$ only. Since we are going to drop $`O(t)`$ and $`O(p^2)`$ terms in the Compton amplitude, we ignore the dependence on $`z^2`$ and treat this matrix element as a function of just two variables $`(pz)`$ and $`(rz)`$. Incorporating the spectral properties of nonforward matrix elements  , we write the plane wave expansion in the form

$`pr/2|𝒪(z|\mathrm{\hspace{0.17em}0})|p+r/2=2(pz){\displaystyle _1^1}𝑑\stackrel{~}{x}{\displaystyle _{1+|\stackrel{~}{x}|}^{1|\stackrel{~}{x}|}}e^{i(kz)}f(\stackrel{~}{x},\alpha )𝑑\alpha +(rz){\displaystyle _1^1}e^{i\alpha (rz)/2}D(\alpha )𝑑\alpha ,`$ (14)
where $`k=\stackrel{~}{x}p+\alpha r/2`$, $`f(\stackrel{~}{x},\alpha )`$ is the double distribution (DD) and $`D(\alpha )`$ is the Polyakov-Weiss (PW) distribution amplitude  absorbing the $`(pz)`$-independent terms. From this parametrization, we can obtain the matrix elements of original uncontracted string operators, (13), including the kinematical twist–3 contributions. We consider first the part coming from the double distribution term in Eq. (14); the contributions from the PW–term will be included separately. In matrix elements, the total derivative turns into the momentum transfer, $`i/X_\sigma r_\sigma =2\xi p_\sigma +\mathrm{\Delta }_\sigma `$. Similarly, we write $`k=(\stackrel{~}{x}+\xi \alpha )p+\alpha \mathrm{\Delta }/2`$. This gives

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_\sigma (z|\mathrm{\hspace{0.17em}0})|p+r/2={\displaystyle _1^1}d\stackrel{~}{x}{\displaystyle _{1+|\stackrel{~}{x}|}^{1|\stackrel{~}{x}|}}d\alpha f(\stackrel{~}{x},\alpha )\{p_\sigma e^{i(\stackrel{~}{x}+\xi \alpha )(pz)i\alpha (\mathrm{\Delta }z)/2}`$ (15)
$`+{\displaystyle \frac{1}{2}}[\mathrm{\Delta }_\sigma (pz)p_\sigma (\mathrm{\Delta }z)]{\displaystyle _0^1}dvve^{iv(\stackrel{~}{x}+\xi \alpha )(pz)iv\alpha (\mathrm{\Delta }z)/2}[\mathrm{sin}(\overline{v}(rz)/2)i\alpha \mathrm{cos}(\overline{v}(rz)/2)]\}.`$ (16)
Skewed distributions. Expanding $`\mathrm{exp}[i\alpha (\mathrm{\Delta }z)/2]=1i\alpha (\mathrm{\Delta }z)/2+\mathrm{}`$ and keeping only terms up to those linear in the transverse momentum***Because of this truncation, the $`\mathrm{\Delta }_\mu \mathrm{\Delta }_\nu `$ terms in the expression for the amplitude $`T_{\mu \nu }`$ will be lost. If needed, they can be kept; see the discussion after Eq. (27). $`\mathrm{\Delta }`$ we get an expression in which the spectral parameter $`\stackrel{~}{x}`$ appears in the exponential factors only in the combination $`x\stackrel{~}{x}+\xi \alpha `$. Thus, we can introduce two skewed parton distributions:

$$\begin{array}{c}\hfill H(x,\xi )\\ \hfill A(x,\xi )\end{array}\}_1^1𝑑\stackrel{~}{x}_{1+|\stackrel{~}{x}|}^{1|\stackrel{~}{x}|}𝑑\alpha \delta (x\stackrel{~}{x}\xi \alpha )f(\stackrel{~}{x},\alpha )\{\begin{array}{c}\hfill 1\\ \hfill \alpha \end{array}$$
(17)
Note that, in our case, the DD $`f(\stackrel{~}{x},\alpha )`$ is even in $`\alpha `$ and odd in $`\stackrel{~}{x}`$. As a result, the functions $`H`$ and $`A`$ satisfy the symmetry relations

$`H(x,\xi )=H(x,\xi ),H(x,\xi )=H(x,\xi ),A(x,\xi )=A(x,\xi ),A(x,\xi )=A(x,\xi ).`$ (18)
Furthermore, because of the antisymmetry of the combination $`\alpha f(\stackrel{~}{x},\alpha )`$ with respect both to $`x`$ and $`\alpha `$ we have

$$_0^1𝑑xA(x,\xi )=\mathrm{\hspace{0.33em}\hspace{0.33em}0}.$$
(19)
Hence, the distribution $`A(x,\xi )`$ cannot be a positive-definite function on $`0x1`$.

Uniting the cosine and sine functions with the overall exponential factor, $`e^{ivx(pz)}`$, one gets $`vx\pm \overline{v}\xi `$ combinations. Using (18), one can arrange that only $`vx+\overline{v}\xi `$ would appear:

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_\sigma (z|\mathrm{\hspace{0.17em}0})|p+r/2=p_\sigma {\displaystyle _1^1}𝑑xe^{ix(pz)}\left[H(x,\xi ){\displaystyle \frac{i(\mathrm{\Delta }z)}{2}}A(x,\xi )\right]`$ (20)
$`+`$ $`{\displaystyle \frac{i}{2}}\left[\mathrm{\Delta }_\sigma (pz)p_\sigma (\mathrm{\Delta }z)\right]{\displaystyle _1^1}𝑑x\left[H(x,\xi )A(x,\xi )\right]{\displaystyle _0^1}𝑑vv\mathrm{cos}[(vx+\overline{v}\xi )(pz)].`$ (21)
In a similar fashion, we get parametrization for the matrix element of the axial string operator (13):

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_{5\sigma }(z|\mathrm{\hspace{0.17em}0})|p+r/2={\displaystyle \frac{i}{2}}ϵ_{\sigma \alpha \beta \gamma }z_\alpha \mathrm{\Delta }_\beta p_\gamma {\displaystyle _1^1}𝑑x\left[H(x,\xi )A(x,\xi )\right]{\displaystyle _0^1}𝑑vv\mathrm{sin}[(vx+\overline{v}\xi )(pz)].`$ (22)
Note that it is expressed in terms of the same skewed distributions $`H(x,\xi )`$ and $`A(x,\xi )`$ which, in turn, are determined by the original double distribution $`f(\stackrel{~}{x},\alpha )`$, see Eq. (17).

## IV  DVCS amplitude for pion target

DD-generated contribution. Substituting the parametrizations (21) and (22) into Eq. (4) and performing the Fourier integral over the separation $`z`$ we obtain the Compton amplitude

$`T_{\mu \nu }`$ $`=`$ $`{\displaystyle \frac{1}{(pq)}}\left[p_\mu q_\nu +q_\mu p_\nu g_{\mu \nu }(pq)+2\xi p_\mu p_\nu +{\displaystyle \frac{\mathrm{\Delta }_\mu }{2}}p_\nu p_\mu {\displaystyle \frac{\mathrm{\Delta }_\nu }{2}}\right]{\displaystyle _1^1}𝑑x{\displaystyle \frac{H(x,\xi )}{x\xi +i0}}`$ (23)
$`+`$ $`{\displaystyle \frac{1}{2(pq)}}{\displaystyle _1^1}𝑑xR(x,\xi ){\displaystyle _0^1}𝑑v{\displaystyle \frac{(q_\mu +3\xi p_\mu )\mathrm{\Delta }_\nu }{\xi +vx+\overline{v}\xi i0}}+{\displaystyle \frac{1}{2(pq)}}{\displaystyle _1^1}𝑑xR(x,\xi ){\displaystyle _0^1}𝑑v{\displaystyle \frac{\mathrm{\Delta }_\mu (q_\nu +\xi p_\nu )}{\xi +vx+\overline{v}\xi +i0}},`$ (24)
where $`R(x,\xi )`$ is a new SPD describing the kinematical twist-3 contributions:

$`R(x,\xi ){\displaystyle \frac{H(x,\xi )}{x}}{\displaystyle \frac{A(x,\xi )}{x}}.`$ (25)
This result for the Compton amplitude contains the same tensor structures as those obtained in Refs.. All three terms in Eq. (24) are individually transverse up to terms of order $`t,p^2`$.

Singularities. The first term is the twist–2 part with the tensor structure corrected exactly as suggested by Guichon and Vanderhaeghen . The integral over $`x`$ exists if $`H(x,\xi )`$ is continuous at $`x=\xi `$, which is the case for SPD’s derived from the DD’s which are less singular than $`1/\stackrel{~}{x}^2`$ for $`\stackrel{~}{x}=0`$ and are continuous otherwise (see ). In particular, continuous SPD’s were obtained in model calculations of SPD’s at a low scale in the instanton vacuum . The second term contributes only to the helicity amplitude for a longitudinally polarized initial photon. The parameter integral over $`v`$ gives the function $`[\mathrm{ln}(x+\xi i0)\mathrm{ln}(2\xi i0)]/(x\xi )`$ which is regular at $`x=\xi `$ and has a logarithmic singularity at $`x=\xi `$. The integral over $`x`$ exists if $`R(x,\xi )`$ is bounded at $`x=\xi `$, which again is the case in the DD-based models described in Ref. . The third term of Eq. (24) corresponds to the transverse polarization of the initial photon. In this case, one faces the integrand $`1/[v(x\xi )+i0]`$ which produces $`dv/v`$ divergence for the $`v`$-integral at the lower limit. One may hope to get a finite result only if the integral

$$I(\xi )_1^1𝑑x\frac{R(x,\xi )}{x\xi +i0}$$
(26)
vanishes. From the definition of the skewed distributions $`H(x,\xi )`$ and $`A(x,\xi )`$ (17) it follows that

$$\frac{A(x,\xi )}{x}=\frac{H(x,\xi )}{\xi }.$$

Hence, one can substitute $`R(x,\xi )`$ by the combination $`H(x,\xi )/x+H(x,\xi )/\xi `$ similar to that used in  within the context of WW approximation (in Refs. , full SPD’s are implied while our $`H(x,\xi )`$ does not include the PW term). By analogy with Ref. , we integrate the $`H(x,\xi )/x`$ term by parts. This gives

$`I(\xi )={\displaystyle \frac{d}{d\xi }}{\displaystyle _1^1}𝑑x{\displaystyle \frac{H(x,\xi )}{x\xi +i0}},`$ (27)
i.e., the $`\xi `$ derivative of the twist-2 contribution. In general, the latter has a nontrivial $`\xi `$-dependent form determined by the shape of SPDs (see, however, the discussion of the PW contribution below). We conclude that the twist-3 part of the tensor amplitude $`T_{\mu \nu }`$ diverges in case of the transverse polarization of the initial photon. A similar observation has been recently made in Ref.. However, it is easy to see that the relevant tensor structure $`\mathrm{\Delta }_\mu (q_\nu +\xi p_\nu )`$ is just a truncated version of the exactly gauge invariant combination $`\mathrm{\Delta }_\mu q_{2}^{}{}_{\nu }{}^{}`$ which has zero projection onto the polarization vector $`ϵ_2^\nu `$ of the final real photon: $`(ϵ_2q_2)=0`$. The structure $`\mathrm{\Delta }_\mu q_{2}^{}{}_{\nu }{}^{}`$ is obtained if one uses the original full form of the DD parametrization (16). It appears from the term with the exponential factor of the argument $`i[v(\stackrel{~}{x}+\xi \alpha )(pz)+\overline{v}\xi (pz)+(v\alpha +\overline{v})(\mathrm{\Delta }z)/2]`$ which is obtained by combining the sine/cosine functions and the exponential factor in the second term of Eq. (16). In the Compton amplitude, it gives rise to a contribution in which the argument of the quark propagator is $`q+v(\stackrel{~}{x}+\xi \alpha )p+\overline{v}\xi p+(v\alpha +\overline{v})\mathrm{\Delta }/2`$. Since $`(\mathrm{\Delta }q)`$, $`(\mathrm{\Delta }p`$) and $`\mathrm{\Delta }^2`$ are negligible, the denominator factors in Eq. (24) remain unchanged. In numerators, representing $`(v\alpha +\overline{v})\mathrm{\Delta }/2`$ as $`[1(1\alpha )v]\mathrm{\Delta }/2`$, we observe that $`\mathrm{\Delta }_\mu (q_\nu +\xi p_\nu )`$ converts into the $`\mathrm{\Delta }_\mu (q_\nu +\xi p_\nu +\mathrm{\Delta }_\nu /2)=\mathrm{\Delta }_\mu q_{2}^{}{}_{\nu }{}^{}`$ term plus a $`v\mathrm{\Delta }_\mu \mathrm{\Delta }_\nu `$ type contribution corresponding to a new SPD built from the $`(1\alpha )^2f(\stackrel{~}{x},\alpha )`$ DD (cf. (17)). Due to the extra $`v`$ factor, the $`v`$-integral for the latter contribution is finite. Hence, for the physical DVCS amplitude, we find no evidence for factorization breaking in the kinematical twist–3 contributions, both in their $`1/\sqrt{q^2}`$ and $`1/q^2`$ terms. It is quite possible that factorization breaks down at the $`1/q^2`$ level, but one needs to analyze $`O(z^2)`$ suppressed terms (i.e., twist–4 contributions) to see if it really happens.

WW-type representation. We can express our results in another form, introducing new skewed distributions related to $`R(x,\xi )`$ $`via`$ an integral transformation similar to that used by Wandzura and Wilczek  WW-type integrals of parton distributions have originally appeared within the $`\xi `$-scaling formalism  . . This allows us to establish the equivalence between our approach and the WW-type approximation proposed in Ref.. Treating the combination $`xv+\overline{v}\xi `$ in (24) as a new variable we define

$`R_W(x,\xi ){\displaystyle _1^1}R(y,\xi )𝑑y{\displaystyle _0^1}\delta (yv+\overline{v}\xi x)𝑑v=\theta (x>\xi ){\displaystyle _x^1}{\displaystyle \frac{R(y,\xi )}{y\xi }}𝑑y\theta (x<\xi ){\displaystyle _1^x}{\displaystyle \frac{R(y,\xi )}{y\xi }}𝑑y.`$ (28)
In terms of this transform, the matrix element of the vector operator (21) can be expressed as

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_\sigma (z|\mathrm{\hspace{0.17em}0})|p+r/2`$ $`=`$ $`{\displaystyle _1^1}dxe^{ix(pz)}\{p_\sigma H(x,\xi ){\displaystyle \frac{i}{2}}p_\sigma (\mathrm{\Delta }z)A(x,\xi )`$ (29)
$`+`$ $`{\displaystyle \frac{1}{4}}(\mathrm{\Delta }_\sigma p_\sigma {\displaystyle \frac{(\mathrm{\Delta }z)}{(pz)}})[R_W(x,\xi )R_W(x,\xi )]\}.`$ (30)
Note that only the odd part of $`R_W(x,\xi )`$ contributes here. In case of the axial operator (22)

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_{5\sigma }(z|\mathrm{\hspace{0.17em}0})|p+r/2`$ $`=`$ $`{\displaystyle \frac{i}{4}}ϵ_{\sigma \alpha \beta \gamma }{\displaystyle \frac{z_\alpha }{(pz)}}\mathrm{\Delta }_\beta p_\gamma {\displaystyle _1^1}𝑑xe^{ix(pz)}\left[R_W(x,\xi )+R_W(x,\xi )\right]`$ (31)
only the even part of $`R_W(x,\xi )`$ appears. The part of the Compton amplitude (24) containing $`R(x,\xi )`$ can be written in terms of this new function as

$`{\displaystyle \frac{1}{2(qp)}}{\displaystyle _1^1}\left[{\displaystyle \frac{\mathrm{\Delta }^\mu (q^\nu +\xi p^\nu )}{x\xi +i0}}+{\displaystyle \frac{(q^\mu +3\xi p^\mu )\mathrm{\Delta }^\nu }{x+\xi i0}}\right]R_W(x,\xi )𝑑x.`$ (32)
The integrals with $`1/(x\pm \xi i0)`$ converge only if the function $`R_W(x,\xi )`$ is continuous for $`x=\pm \xi `$. According to Eq. (28), $`R_W(x,\xi )`$ is given by the integral of $`R(y,\xi )/(y\xi )`$ from $`x`$ to 1 if $`x>\xi `$ and from $`x`$ to $`1`$ if $`x<\xi `$. Evidently, $`x=\xi `$ is not a special point in the integral transformation (28), hence the function $`R_W(x,\xi )`$ is continuous at $`x=\xi `$. However, it is extremely unlikely that the limiting values approached by $`R_W(x,\xi )`$ for $`x=\xi `$ from below and from above do coincide. Indeed, the difference of the two limits can be written as the principal value integral (compare with )

$`R_W(\xi +0,\xi )R_W(\xi 0,\xi )=\mathrm{P}{\displaystyle _1^1}{\displaystyle \frac{R(y,\xi )}{y\xi }}𝑑y,`$ (33)
which can be converted into the $`\xi `$-derivative of the real part of the twist–2 contribution. This means that the singularity, which we observed as a straight divergence of the $`dv/v`$ integral, in this approach appears due to an unavoidable discontinuity of the $`R_W(x,\xi )`$ transform at $`x=\xi `$.

Contribution from the PW–term. The contribution of the PW term to the vector operator

$`{\displaystyle \frac{1}{2}}pr/2|𝒪_\sigma (0|z)|p+r/2_{\text{PW–Term}}={\displaystyle \frac{r_\sigma }{2}}{\displaystyle _1^1}𝑑\alpha e^{i\alpha (rz)/2}D(\alpha )`$ (34)
has a simple structure corresponding to a parton picture in which the partons carry the fractions $`(1\pm \alpha )/2`$ of the momentum transfer $`r`$. Since only one momentum $`r`$ is involved, this term can contribute only to the totally symmetric part of the vector string operator: it “decouples” in the reduction relations (11). In particular, the PW term does not contribute to the second contribution in Eq. (13) which is generated by decomposition of the axial string operator: both derivatives, with respect to $`X`$ and $`z`$, give rise to the momentum transfer $`r`$, whence the contraction with the $`ϵ`$–tensor in (13) gives zero. Thus, the PW-contribution should be transverse by itself. Indeed, a straightforward calculation gives

$$T_{\mu \nu }|_{\text{PW}}=\frac{1}{(rq)}\left[r_\mu q_\nu +q_\mu r_\nu g_{\mu \nu }(rq)+r_\mu r_\nu \right]_1^1\frac{D(\alpha )}{\alpha 1}𝑑\alpha ,$$
(35)
which evidently satisfies $`q_\mu T_{\mu \nu }|_{\text{PW}}=0,r_\mu T_{\mu \nu }|_{\text{PW}}=0.`$ Hence, this term can be treated as a separate contribution.

Alternatively, one may include it into the basic SPD $`H(x,\xi )`$ and all SPD’s derived from $`H(x,\xi )`$. Specifically, for $`\xi >0`$, the PW contribution to $`H(x,\xi )`$ is $`D(x/\xi )\theta (|x|\xi )`$ ; it contributes $`(\xi x)D^{}(x/\xi )\theta (|x|\xi )/\xi ^2`$ \[where $`D^{}(\alpha )(d/d\alpha )D(\alpha )]`$ to $`R(x,\xi )`$; furthermore, the PW contribution to $`R_W(x,\xi )`$ is $`D(x/\xi )\theta (|x|\xi )/\xi `$. Inserting these functions into Eqs. (24) and (32) one rederives Eq. (35). One can also observe that the PW term gives zero contribution into $`I(\xi )`$, Eq. (26).

## V Conclusions

In this paper, we have studied the DVCS amplitude $`T_{\mu \nu }`$ making use of the light–cone expansion in terms of QCD string operators in the coordinate space. We have demonstrated that transversality of the light–cone expansion can be maintained by including a minimal set of “kinematical” twist–3 operators, which appear as total derivatives of twist–2 operators. Incorporating the formalism of double distributions, we established that the kinematical twist–3 contributions are described by a skewed parton distribution $`R_W(x,\xi )`$ which can be derived from the basic twist–2 double distribution $`f(\stackrel{~}{x},\alpha )`$. The new SPD $`R_W(x,\xi )`$ has the structure of a generalized WW transform. Due to discontinuities of $`R_W(x,\xi )`$ at $`x=\xi `$, the factorization for the tensor amplitude $`T_{\mu \nu }`$ breaks down at the twist-3 level for the part corresponding to the transverse virtual photon. However, this term gives zero projection onto the polarization vector of the final real photon.

## VI Acknowledgements

We are grateful to I.I. Balitsky, J. Blümlein, V.M. Braun, B. Geyer, N. Kivel, M. Lazar, I.V. Musatov, M.V. Polyakov and O.V. Teryaev for many stimulating discussions. We thank M. Diehl for useful comments. This work was supported by the US Department of Energy under contract DE-AC05-84ER40150, by the Deutsche Forschungsgemeinschaft (DFG), by the German Ministry of Education and Research (BMBF), and by COSY, Jülich. A.R. expresses his deep gratitude to K. Goeke for warm hospitality in Bochum and support.


text:
# Effective attraction between like-charged colloids in a 2D plasma

## I Introduction

Recent experiments show convincingly the existence of attractions between like-charged colloids immersed in ionic solution, in particular in the vicinity of a glass wall or when the colloids are confined between glass walls . This remarkable counter-intuitive phenomenon is inconsistent with the well established theory of Derjaguin, Landau, Verwey and Overbeek (DLVO) and has generated various theoretical interpretations.

It was shown  that counterion correlation forces, because of the long range of the electrostatic potentials, are usually not pairwise additive. Many-body effects  and Coulomb depletion forces  have been invoked to explain attractive forces between like-charged macroions. An exact demonstration  proved that the non-linear Poisson-Boltzmann (PB) mean-field equation can not give attractions in the case of Dirichlet boundary conditions. The argument was extended to a broader class of models  and boundary conditions . Recently, there was a significant new attempt based on non-equilibrium hydrodynamic effects .

In this paper, we solve a 2D one-component plasma model exactly at a certain coupling constant following the method of Jancovici . We have found a semi-analytic method which extends the exact solution to the case in which colloidal particles are present. We treat the colloids in the plasma as empty regions (voids) from which the ions are excluded. The colloid may have fixed charge at its center. The analytic results obtained for small-sized colloids show that while the ion-ion repulsion, which goes like $`\mathrm{e}^{R^2}`$ ($`R`$ is the separation distance) in this particular 2D model, is strong at short distances, the ion-void attraction, which goes like $`R^2\mathrm{e}^{R^2}`$, becomes dominant at large distances. The latter property is not normally found in conventional mean field approximations, and must come from strong correlation effects. For large-sized colloids, we also have numerical results which qualitatively agree with the analytic ones. The present state of understanding of how to treat reliably fluctuations beyond mean field theory is not very good despite the fact that this is an old problem . Thus, our exact results may serve as a testing ground for more sophisticated approximation schemes. In addition, since this is a 2D model, the results are of interest in connection with the problem of interactions among linear polyelectrolyte molecules .

In the sections that follow, we shall first review the plasma model in section II. Then in section III and IV we calculate the ion-ion and void-void interactions in the plasma system respectively. Section V puts the above results together by treating the colloid as an empty region plus a fixed ion charge at the center. We have relegated the technical details as much as possible to the appendices.

## II The Plasma System

We start by reviewing the exact solution for the one-component 2D plasma system . Our approach is to map the classical plasma system into a quantum Hall system. Such an analogy allows one to describe the interacting ions in terms of non-interacting electrons occupying certain quantum obitals, and subsequently use already developed tools. The advantage of this approach will become more clear later when we extend the same analogy to more complicated systems with fixed ion charges or/and empty regions present.

### A The Quantum Orbital Language

The one-component 2D plasma system is composed of $`n`$ simple ions embedded in a disk of uniform background charge. Each ion carries a charge of $`1`$; the disk has a radius $`r`$ and an areal charge density of $`\rho _0`$. Throughout this paper, we shall express all lengths in the unit of $`\mathrm{}_0=\frac{1}{\sqrt{\pi \rho _0}}`$. Thus, a unit circle contains unit amount of background charge, and the charge neutrality condition can be written simply as

$$n=r^2.$$
(1)
In two dimensions, it is convenient to use a complex number $`z=x+𝐢y`$ (in units of $`\mathrm{}_0`$) to represent the location of $`(x,y)`$. Using $`z_i`$ to denote the position of the $`i`$th ion, we can solve the 2D Poisson equation and find the total Coulomb potential energy of the system:

$$V=\frac{1}{2}\underset{i=1}{\overset{n}{}}|z_i|^2\underset{i<j}{\overset{n}{}}\mathrm{ln}|z_iz_j|.$$
(2)
Here the first term describes the attractions between the ions and the background charge, and the second term represents the mutual repulsions among the ions themselves.

At inverse temperature $`\beta =2`$, the Boltzmann factor reads:

$$\mathrm{e}^{\beta V}=\mathrm{e}^{_{i=1}^n|z_i|^2}\left|\underset{i<j}{\overset{n}{}}(z_iz_j)\right|^2.$$
(3)
Written in the above form, the Boltzmann factor can be identified as the square of the determinant of a matrix $`𝐌`$

$$\mathrm{e}^{\beta V}=|det𝐌|^2,$$
(4)
with the matrix elements being the lowest Landau level (LLL) wavefunctions studied in connection with the quantum Hall effect :

$$M_{ij}=\psi _i(z_j)\frac{z_j^{i1}}{\sqrt{\pi \mathrm{\Gamma }(i)}}\mathrm{e}^{|z_j|^2/2}.$$
(5)
The above wavefunction describes an electron occupying the $`i`$th angular momentum orbital in the LLL. In this quantum language, $`det𝐌`$ is a Slater determinant and represents the Laughlin wavefunction at filling fact $`\nu =1`$, describing a fully filled Landau level. Thus, the complicated Coulomb interactions are replaced by a simple Slater determinant. The determinant representation of the Boltzmann factor was first noted by Jancovici  prior to the discovery of the quantum Hall effect. We recapitulate (and then extend) his argument in more modern language.

### B The Partition Function

Armed with this powerful quantum analogy, we proceed to calculate the partition function of the system, by integrating out all the ionic degrees of freedom:

$$Z=\underset{i=1}{\overset{n}{}}dz_i\mathrm{e}^{\beta V}.$$
(6)
This is nothing more than the norm of the Laughlin wavefunction, and hence can be computed exactly. Substituting eq.(4) in eq.(6) and expanding the determinant, we have:

$`Z`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{n}{}}}{\displaystyle dz_i|det𝐌|^2}`$ (7)
$`=`$ $`{\displaystyle \underset{i=1}{\overset{n}{}}}{\displaystyle dz_i\underset{\{𝒫\}}{}\underset{\{𝒬\}}{}(1)^{𝒫+𝒬}\underset{i=1}{\overset{n}{}}\overline{\psi }_{P_i}(z_i)\psi _{Q_i}(z_i)},`$ (8)
where $`𝒫`$ and $`𝒬`$ are permutations of $`\{1:n\}`$ \[$`1:n`$ is a short hand notation for $`1,2,\mathrm{},n`$\] that come from the expansions of $`det\overline{𝐌}`$ and $`det𝐌`$ respectively. In eq.(8), we have both summations and integrations which, respectively represent the Coulomb interactions and the thermal averages. It is useful at this point to postpone treating the Coulomb interactions until after the thermal averages are done. Thus, we switch the order of the summations and integrations in eq.(8) so that:

$`Z`$ $`=`$ $`{\displaystyle \underset{\{𝒫\}}{}}{\displaystyle \underset{\{𝒬\}}{}}(1)^{𝒫+𝒬}{\displaystyle \underset{i=1}{\overset{n}{}}}{\displaystyle dz_i\overline{\psi }_{P_i}(z_i)\psi _{Q_i}(z_i)}`$ (9)
$``$ $`{\displaystyle \underset{\{𝒫\}}{}}{\displaystyle \underset{\{𝒬\}}{}}(1)^{𝒫+𝒬}{\displaystyle \underset{i=1}{\overset{n}{}}}P_i|Q_i,`$ (10)
where $`i|j`$ is the Dirac notation for the inner product of the $`i`$th and $`j`$th LLL orbitals. Notice that the summations over $`\{𝒬\}`$ can be identified as a matrix determinant:

$$\underset{\{𝒬\}}{}(1)^{𝒫+𝒬}\underset{i=1}{\overset{n}{}}P_i|Q_i=\underset{\{\}}{}(1)^{}\underset{i=1}{\overset{n}{}}i|R_i=det𝐎,$$
(11)
where $`𝒬𝒫^1`$ is the composite permutation of $`𝒬`$ and $`𝒫^1`$, and the overlap matrix $`𝐎`$ defined by

$$O_{ij}=i|j$$
(12)
is a matrix composed of LLL wavefunction inner products. Here, just as previously, the complicated Coulomb interactions are again replaced by a simple matrix determinant. Finally, we substitute eq.(11) into (10) to obtain:

$$Z=\underset{\{𝒫\}}{}det𝐎=n!det𝐎.$$
(13)
The immaterial factor of $`n!`$ will be dropped henceforth. Through the quantum analogy, the statistical problem is reduced to evaluating the overlap matrix $`𝐎`$ and its determinant. This is particularly simple in the present case, because of orthonormality we know that matrix $`𝐎`$ is an identity matrix.

### C The Plasma Edge

Near the plasma edge, physical properties are quite different from those in the bulk. This motivates us to study the boundary effects separately. For simplicity, we are going to be mainly interested in a soft boundary condition, which restricts the range of the ions’ angular momenta rather than their positions. The ions can go anywhere in the complex plane, nevertheless, they are most likely to be found inside a disk with radius $`r=\sqrt{n}`$ because their angular momenta are bounded by $`n`$ \[Recall that the $`l`$th LLL state has angular momentum of $`(l1)`$ and peaks in a shell $`\sqrt{l1}<r<\sqrt{l}`$\]. The advantage of using this type of boundary condition is apparent: we can continue to use the same wavefunctions to describe the ions near the edge. What is implicitly assumed here is that the uniform background charge does not cease at the plasma edge; it extends to infinity. We also note that in the thermodynamic limit ($`r\mathrm{}`$), the circular plasma edge becomes locally flat.

One may also consider putting some ‘surface’ charge on the plasma boundary. But, as indicated below, this effect can be absorbed in shifting the position of the boundary. We assume the surface charge, if any, does not fluctuate. Depending on its sign, the boundary charge either increases or decreases the number of available angular momentum channels. This is more or less equivalent to varying the radius $`r`$ of the plasma disk in a system with neutral boundary. Hence, without loss of generality, we may assume boundaries are neutral.

Another type of plasma system we want to address in this paper is a plasma strip with two parallel edges separated by $`h`$. This can be realized in the thermodynamic limit by a system of ions with angular momenta ranging from $`(r\frac{h}{2})^2`$ to $`(r+\frac{h}{2})^2`$. The double edged system thus differs from single edged system only in the momenta range. To avoid redundancy, we shall use the single edged system as example in our calculations. By changing the angular momentum limits, the results can be easily transferred to the double edged systems.

## III Fixed Ion Charges in a Plasma

In this section, we shall apply the same method to a slightly more complicated plasma system, with one or two ion charge(s) fixed at certain position(s). The resulting partition function is conventionally called the one- or two-body density function. The physical significance is two-fold. First, we can learn from these density functions about how the ions are distributed and how two ions are correlated. More importantly, the fixed ion charges can be viewed as the limiting case of small-sized charged colloids. Thus, the results obtained in this section will provide us with many useful insights to the colloidal interactions.

### A The One-body Density Function

First we consider the one-body density function. We fix the position $`z_1=w`$ and integrate out the remaining ionic degrees of freedom

$`Z^{(1)}`$ $`=`$ $`{\displaystyle \underset{i=2}{\overset{n}{}}}{\displaystyle dz_i\mathrm{e}^{\beta V}}|_{z_1=w}`$ (14)
$`=`$ $`{\displaystyle \underset{i=2}{\overset{n}{}}}{\displaystyle dz_i|det𝐌|^2}|_{z_1=w}`$ (15)
$`=`$ $`{\displaystyle \underset{\{𝒫\}}{}}{\displaystyle \underset{\{𝒬\}}{}}(1)^{𝒫+𝒬}\overline{\psi }_{P_1}(w)\psi _{Q_1}(w){\displaystyle \underset{i=2}{\overset{n}{}}}P_i|Q_i.`$ (16)
The summations over $`P_1`$ and $`\{𝒬\}`$ then reduce to a matrix determinant:

$$\underset{P_1}{}\underset{\{𝒬\}}{}(1)^{𝒫+𝒬}\overline{\psi }_{P_1}(w)\psi _{Q_1}(w)\underset{i=2}{\overset{n}{}}P_i|Q_i=det𝐎^{(1)},$$
(17)
with the $`(n+1)`$ by $`(n+1)`$ matrix $`𝐎^{(1)}`$ defined by:

$$𝐎^{(1)}\left[\begin{array}{cc}𝐎& \mathrm{\Psi }^{}(w)\\ \mathrm{\Psi }(w)& 0\end{array}\right].$$
(18)
Here, $`𝐎`$ (not to be confused with $`𝐎^{(1)}`$) is the $`n\times n`$ overlap matrix defined in (12); $`\mathrm{\Psi }(\psi _1,\psi _2,\mathrm{},\psi _n)`$ is a row vector of LLL wavefunctions representing the fixed ion charge, and $`\mathrm{\Psi }^{}`$ is its conjugate transpose. Finally, the remaining summations over $`\{P_{i,(i=2:n)}\}`$ in eq.(16) just yield a trivial prefactor:

$$Z^{(1)}=\underset{\{P_{i(i=2:n)}\}}{}det𝐎^{(1)}=(n1)!det𝐎^{(1)}.$$
(19)
To evaluate the determinant of matrix $`𝐎^{(1)}`$, we apply Schur’s theorem (discussed in detail in appendix B) and find:

$`Z^{(1)}`$ $`=`$ $`det𝐎^{(1)}=det𝐎det\left[\mathrm{\Psi }(w)𝐎^1\mathrm{\Psi }^{}(w)\right]_{1\times 1}`$ (20)
$`=`$ $`{\displaystyle \underset{k=1}{\overset{n}{}}}\psi _k(w)\overline{\psi }_k(w),`$ (21)
where we have used the fact that the overlap matrix $`𝐎`$ is an identity matrix. Using eq.(A5), we convert the above series into translation coefficients (TCs) defined in eq.(A4):

$$Z^{(1)}=\underset{k=1}{\overset{n}{}}\psi _k(w)\overline{\psi }_k(w)=\frac{1}{\pi }\underset{k=1}{\overset{n}{}}|T_{k1}(w)|^2.$$
(22)
Now that the partition function is written in the form of a TC sum, we may use the TC sum-rules listed in appendix A 2. We consider the following two cases.

If the fixed ion is in the bulk of the plasma ($`wr=\sqrt{n}`$), we may take the limit $`n\mathrm{}`$ first and use TC sum-rule (A20) to obtain:

$$Z^{(1)}=\frac{1}{\pi }.$$
(23)
This tells that deep in the bulk, the ion density is constant and equal to the background charge density everywhere.

Near the plasma edge, we assume the boundary is located at $`(r+\frac{h}{2})`$ and the fixed ion is located at $`w=r\mathrm{e}^{𝐢\theta }`$, so that the ion is $`\frac{h}{2}`$ away from the boundary. In this case, we can use TC sum-rule (A30) in eq.(22) to obtain:

$$Z^{(1)s}=\frac{1}{2\pi }\left[1+\mathrm{\Phi }\left(\frac{h}{\sqrt{2}}\right)\right].$$
(24)
Here, the superscript ‘s’ means single edged system (in contrast to ‘d’, which means double edged system) and $`\mathrm{\Phi }(x)`$ is the standard error function. As shown in Fig. 1, the ion density indeed drops from its bulk value to zero upon crossing the soft plasma boundary (at $`h=0`$).

### B The Two-body Density Function

In calculating the one-body density function, we find that a fixed ion charge amounts to an additional row and column to the overlap matrix. This result can be generalized for systems with multiple fixed ion charges. In particular, we consider the two-body density function with $`z_1=w`$ and $`z_2=\overline{w}`$ fixed. The result is:

$$Z^{(2)}=det𝐎^{(2)},$$
(25)
where $`𝐎^{(2)}`$ is an $`(n+2)`$ by $`(n+2)`$ matrix:

$$𝐎^{(2)}\left[\begin{array}{ccc}𝐎& \mathrm{\Psi }^{}(w)& \mathrm{\Psi }^{}(\overline{w})\\ \mathrm{\Psi }(w)& 0& 0\\ \mathrm{\Psi }(\overline{w})& 0& 0\end{array}\right].$$
(26)
Similarly using Schur’s theorem and eq.(A5), we can write the correlation function in terms of TC sums:

$`Z^{(2)}`$ $`=`$ $`det𝐎det\left(\left[\begin{array}{c}\mathrm{\Psi }(w)\\ \mathrm{\Psi }(\overline{w})\end{array}\right]𝐎^1\left[\begin{array}{cc}\mathrm{\Psi }^{}(w)& \mathrm{\Psi }^{}(\overline{w})\end{array}\right]\right)_{2\times 2}`$ (30)
$`=`$ $`{\displaystyle \frac{1}{\pi ^2}}\left[{\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}(w)|^2\right]^2{\displaystyle \frac{1}{\pi ^2}}\left|{\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}^2(w)\right|^2.`$ (31)
In the bulk of the plasma ($`wr`$), we assume $`w=𝐢\frac{R}{2}`$ so that the two fixed ions are $`R`$ apart. Plugging TC sum-rules (A20) and (A21) into eq.(31), and using eq.(A15) for $`T_{11}(𝐢R)`$, we obtain:

$$Z^{(2)}=\frac{1}{\pi ^2}\left[1|T_{11}(𝐢R)|^2\right]=\frac{1}{\pi ^2}\left(1\mathrm{e}^{R^2}\right).$$
(32)
The above result tells that the Coulomb repulsion between the two ion charges is heavily screened by the other ions and becomes short ranged. Its decay has a Gaussian form, which is much faster than the exponential form that DLVO theory predicts.

Near the plasma edge, we assume the boundary is located at $`(r+\frac{h}{2})`$ and the fixed ions are located at $`w=r\mathrm{e}^{𝐢\theta }`$ and $`\overline{w}`$, where $`\theta \frac{R}{2r}`$ so that the fixed ions are $`\frac{h}{2}`$ away from the edge and $`R`$ apart from each other. Using TC sum-rules (A30) and (A33) in eq.(31), we obtain (for asymptotically large $`R`$):

$$Z^{(2)s}=\frac{1}{4\pi ^2}\left\{\left[1+\mathrm{\Phi }\left(\frac{h}{\sqrt{2}}\right)\right]^2\frac{2\mathrm{e}^{h^2}}{\pi R^2}\right\}.$$
(33)
The above result shows that the ion-ion correlation function has a power law decay near the plasma edge. This is fundamentally different from the bulk behavior, and is due to the dipole moments induced by the sharp cutoff in angular momentum.

Finally, in the double edged system where the plasma edges are located at $`(r\pm \frac{h}{2})`$, we may use TC sum-rules (A37) and (A40) in eq.(31) to obtain (for asymptotically large $`R`$):

$$Z^{(2)d}=\frac{1}{\pi ^2}\left[\mathrm{\Phi }^2\left(\frac{h}{\sqrt{2}}\right)\frac{\mathrm{e}^{h^2}}{2\pi R^2}\mathrm{sin}^2(hR)\right].$$
(34)
Here, we find the correlation is not only long ranged, but oscillating with a period of $`2\pi h^1`$ as well.

### C Divalent Ions

A divalent ion can be viewed as composed of two simple ions occupying the same place. Here we want to calculate the two-body density function for a pair of divalent ions embedded in the monovalent plasma. Let us assume $`z_1=z_2=w`$ and $`z_3=z_4=\overline{w}`$ are fixed, where $`w=𝐢\frac{R}{2}`$ so that the two divalent ions are $`R`$ apart. Naively one may follow eq.(26) and write down an $`(n+4)`$ by $`(n+4)`$ matrix:

$$\stackrel{~}{𝐎}^{(4)}\left[\begin{array}{ccccc}𝐎& \mathrm{\Psi }^{}(z_1)& \mathrm{\Psi }^{}(z_2)& \mathrm{\Psi }^{}(z_3)& \mathrm{\Psi }^{}(z_4)\\ \mathrm{\Psi }(z_1)& 0& 0& 0& 0\\ \mathrm{\Psi }(z_2)& 0& 0& 0& 0\\ \mathrm{\Psi }(z_3)& 0& 0& 0& 0\\ \mathrm{\Psi }(z_4)& 0& 0& 0& 0\end{array}\right].$$
(35)
However, the matrix $`\stackrel{~}{𝐎}^{(4)}`$ defined above is singular because $`z_1=z_2`$ (or $`z_3=z_4`$). This singularity actually originates from the Boltzman factor in eq.(3), where we incorrectly include the self interaction between $`z_1`$ and $`z_2`$ (or $`z_3`$ and $`z_4`$). To get around this problem, we must extract the singular factors from the matrix determinant. More precisely, we want to consider the following limits:

$$det𝐎^{(4)}\underset{z_2z_1}{lim}\underset{z_4z_3}{lim}\frac{det\stackrel{~}{𝐎}^{(4)}(z_1,z_2,z_3,z_4)}{|z_2z_1|^2|z_4z_3|^2}.$$
(36)
To construct matrix $`𝐎^{(4)}`$ from $`\stackrel{~}{𝐎}^{(4)}`$, we subtract $`\mathrm{\Psi }(z_1)`$ in the $`(n+1)`$th row from $`\mathrm{\Psi }(z_2)`$ in the $`(n+2)`$th row in eq.(35). The resulting $`(n+2)`$th row then yields a factor of $`(z_2z_1)`$:

$$\underset{z_2z_1}{lim}[\mathrm{\Psi }(z_2)\mathrm{\Psi }(z_1)]=(z_2z_1)\mathrm{\Psi }^{}(z_1),$$
(37)
where $`\mathrm{\Psi }^{}(\psi _1^{},\psi _2^{},\mathrm{},\psi _n^{})`$, and $`\psi ^{}`$ is the first derivative of $`\psi `$. More precisely, following the usual procedure for the Hilbert space of analytic functions , we take the derivative to act only on the analytic part of $`\psi `$ and not on the Gaussian factor. The latter represents the interaction with the background charge and does not possess the self-interaction singularity.

The factors of $`(\overline{z}_2\overline{z}_1)`$ and $`|z_4z_3|^2`$ can be extracted via similar procedures and we are left with a matrix that looks like the following:

$$𝐎^{(4)}\left[\begin{array}{ccccc}𝐎& \mathrm{\Psi }^{}(w)& \mathrm{\Psi }_{}^{}{}_{}{}^{}(w)& \mathrm{\Psi }^{}(\overline{w})& \mathrm{\Psi }_{}^{}{}_{}{}^{}(\overline{w})\\ \mathrm{\Psi }(w)& 0& 0& 0& 0\\ \mathrm{\Psi }^{}(w)& 0& 0& 0& 0\\ \mathrm{\Psi }(\overline{w})& 0& 0& 0& 0\\ \mathrm{\Psi }^{}(\overline{w})& 0& 0& 0& 0\end{array}\right].$$
(38)
Applying Schur’s theorem and using various TC sum-rules listed in appendix A 2, one finds:

$`Z^{(4)}`$ $`=`$ $`det𝐎^{(4)}`$ (39)
$`=`$ $`{\displaystyle \frac{1}{\pi ^4}}det\left[\begin{array}{cccc}1& 0& \mathrm{e}^{R^2/2}& 0\\ 0& 1& R\mathrm{e}^{R^2/2}& \mathrm{e}^{R^2/2}\\ \mathrm{e}^{R^2/2}& R\mathrm{e}^{R^2/2}& 1& R\\ 0& \mathrm{e}^{R^2/2}& R& 1+R^2\end{array}\right]`$ (44)
$`=`$ $`{\displaystyle \frac{1}{\pi ^4}}\left[\left(1\mathrm{e}^{R^2}\right)^2R^4\mathrm{e}^{R^2}\right].`$ (45)
To compare the correlation functions of divalent ions and simple ions (in the bulk case), we plot eq.(32) and (45) in Fig. 2. Clearly, the divalent ions has a larger ‘exchange hole’ near its origin.

## IV Void Colloids in a Plasma

In this section, we introduce two empty colloidal voids into the plasma system. We take the colloidal particle to be a circular disk with finite radius $`a`$, inside which ions are excluded. The presence of the voids clearly destroys the system’s azimuthal symmetry. As a result, different angular momentum orbitals intersect, and the overlap matrix is generally not diagonal. It is the off-diagonal matrix elements that make the colloidal interactions non-trivial.

We shall fix the two voids at $`w`$ and $`\overline{w}`$, respectively. Following eq.(13), we can write down the partition function of the system as the determinant of the overlap matrix:

$$Z=det𝐎.$$
(46)
Because the regions occupied by colloids are not accessible to the ions, their contribution to the overlap matrix should be excluded. Thus, in contrast to eq.(12), here we have:

$$O_{ij}=i|ji|j_wi|j_{\overline{w}},$$
(47)
where $`i|j_w`$ means a partial inner product integrated over a circle with radius $`a`$ centered at $`w`$. According to eq.(A9), eq.(47) can be explicitly evaluated as :

$$O_{ij}=\delta _{ij}[\underset{l=1}{\overset{\mathrm{}}{}}\gamma _l\overline{T}_{il}(w)T_{jl}(w)+\mathrm{c}.\mathrm{c}.],$$
(48)
where $`\gamma _l`$ is the incomplete gamma function defined in eq.(A3).

### A Analytic Results

The off-diagonal matrix elements make the determinant in eq.(46) difficult to analyse. However, if the colloidal size $`a`$ is small, we can treat them perturbatively. From eq.(A3) we know $`\gamma _iO(a^{2i})`$ when $`a0`$. Thus, to the lowest order in $`a`$, we can truncate the infinite summation in eq.(48) at $`l=1`$:

$$O_{ij}\delta _{ij}\gamma _1\left[\overline{T}_{i1}(w)T_{j1}(w)+T_{i1}(w)\overline{T}_{j1}(w)\right].$$
(49)
If we define a $`2\times n`$ matrix $`𝐉`$ as

$$𝐉\left[\begin{array}{cccc}T_{11}(w)& T_{21}(w)& \mathrm{}& T_{n1}(w)\\ \overline{T}_{11}(w)& \overline{T}_{21}(w)& \mathrm{}& \overline{T}_{n1}(w)\end{array}\right],$$
(50)
eq.(49) may be rewritten in a matrix form:

$$𝐎=1𝐉^{}\gamma _1𝐉.$$
(51)
According to the corollary of Schur’s theorem (see appendix B, eq.(B3)), we have:

$$det\left(1𝐉^{}\gamma _1𝐉\right)_{n\times n}=det\left(1𝐉\gamma _1𝐉^{}\right)_{2\times 2}.$$
(52)
The determinant of a $`2\times 2`$ matrix is easy to calculate and the result is (dropping terms with order higher than $`\gamma _1`$):

$$Z=det𝐎12\left(\underset{k=1}{\overset{n}{}}|T_{k1}(w)|^2\right)\gamma _1=12\gamma _1,$$
(53)
where we have used the TC sum-rule (A20). Notice that the above (first order) result is $`R`$-independent. This is no accident, because in order for the colloids to know the separation distance $`R`$, the Green’s function needs to be integrated over both colloidal regions, which yeilds a term proportional to $`a^4`$, appearing in at least the second order in perturbation theory.

In the second order approximation, we must keep terms proportional to $`\gamma _1`$, $`\gamma _1^2`$ and $`\gamma _2`$. The procedure is quite similar to the above, and the result is:

$`Z`$ $``$ $`12\left({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}(w)|^2\right)\gamma _12\left({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k2}(w)|^2\right)\gamma _2`$ (55)
$`+\left[\left({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}(w)|^2\right)^2\left|{\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}^2(w)\right|^2\right]\gamma _1^2.`$
Compared to eq.(22) and (31), we find that the first and second order terms are respectively proportional to the one- and two-body density functions. Hence the colloidal void-void interaction does not differ very much from the ion-ion interaction we studied earlier in Section III.

In the plasma bulk, we assume $`w=𝐢\frac{R}{2}`$ and use TC sum-rules (A20) and (A21) in eq.(55) to obtain

$$Z12\gamma _12\gamma _2+\left(1\mathrm{e}^{R^2}\right)\gamma _1^2;$$
(56)
near the plasma edge(s), we assume $`w=r\mathrm{e}^{𝐢\theta }`$ ($`\theta =\frac{R}{2r}`$) and express eq.(55) in terms of TC sum-rules results:

$$Z^{s,d}12A_{11}^{s,d}\gamma _12A_{22}^{s,d}\gamma _2+\left[(A_{11}^{s,d})^2|B_{11}^{s,d}|^2\right]\gamma _1^2,$$
(57)
where the explicit forms for $`A^{s,d}`$ and $`B^{s,d}`$ for asymptotically large $`R`$ can be found in eq.(A30), (A33), (A37) and (A40).

### B Numerical Results

For large-sized colloids, the analytic results based on the small $`a`$-expansion may not be appropriate. We resort to numerical methods and use LU factorization with partial pivoting to evaluate the matrix determinants. We find that (essentially) exact results can be obtained from relatively small matrix sizes.

For colloidal voids immersed in the bulk of the plasma, our numerical results are presented in Fig. 3. For a variety of sizes we studied, no attraction is found, confirming our analytic predictions.

Near the plasma edge(s), we find it more convenient to use square shaped colloids in numerical calculation. We assume the squares are of size $`2a`$, centered at $`r\mathrm{e}^{𝐢\theta }`$ and $`r\mathrm{e}^{𝐢\theta }`$ respectively. \[The plasma edge(s) are located at $`(r\pm \frac{h}{2})`$.\] Defining $`\delta \frac{a}{r}`$, we can calculate the overlap matrix in eq.(47) as following:

$$O_{ij}=\delta _{ij}[_{ra}^{r+a}\mathrm{d}\rho _{\theta \delta }^{\theta +\delta }\mathrm{d}\varphi \frac{\rho ^{i+j1}\mathrm{e}^{\rho ^2}\mathrm{e}^{𝐢(ji)\varphi }}{\pi \sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(j)}}+\mathrm{c}.\mathrm{c}.].$$
(58)
In the limit $`ijr^2\mathrm{}`$, the above integral can be evaluated asymptotically. The result is

$$O_{ij}\delta _{ij}C(\theta ,\delta )\frac{\mathrm{\Gamma }(\frac{i+j}{2})}{\sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(j)}}\left[\mathrm{\Phi }\left(r_+\right)\mathrm{\Phi }\left(r_{}\right)\right],$$
(59)
where, $`r_\pm \frac{(r\pm a)^2\frac{i+j}{2}+1}{\sqrt{i+j2}}`$; $`C(\theta ,\delta )\frac{\mathrm{cos}[(ij)\theta ]\mathrm{sin}[(ij)\delta ]}{(ij)\pi }`$ ($`C(\theta ,\delta )`$ should be interpreted as $`\frac{\delta }{\pi }`$ in case of $`i=j`$). In practice, we keep increasing $`r`$ until the numerical result becomes independent of $`r`$, which is found to occur approximately at $`r=60`$. The results shown in Fig. 4 (for single edged system) and Fig. 5 (for double edged system) are computed using $`r=100`$. In both figures, when the colloidal voids are far away from the edge(s) ($`h=6.0`$), the curves are similar to those in Fig. 3, namely, the repulsive interaction has a Gaussian decay. As the colloids become closer to the plasma edge(s) (smaller $`h`$), the repulsion starts to develop a long tail, in agreement with our analytic results.

## V Charged Colloids in a Plasma

We now consider charged colloids in the plasma. In the simplest case, the charged colloid is taken to be an empty disk plus a fixed ion charge at the center. The ion-ion interaction and void-void interaction have been studied in Section III and IV respectively; both are repulsive. However, as we shall see in this section, the ion-void interaction is usually attractive, and has a longer range. The competition between the repulsion and attraction results in a richer behavior for the charged colloidal interactions.

Following eq.(25) and (26), we have the partition function as follows:

$`Z^{(2)}`$ $`=`$ $`det𝐎^{(2)}`$ (60)
$`𝐎^{(2)}`$ $``$ $`\left[\begin{array}{ccc}𝐎& \mathrm{\Psi }^{}(w)& \mathrm{\Psi }^{}(\overline{w})\\ \mathrm{\Psi }(w)& 0& 0\\ \mathrm{\Psi }(\overline{w})& 0& 0\end{array}\right],`$ (64)
where, $`\mathrm{\Psi }(\psi _1,\psi _2,\mathrm{},\psi _n)`$ is a row vector representing the fixed ion charge at the colloidal center; the $`n\times n`$ overlap matrix $`𝐎`$, according to eq.(48), is non-diagonal:

$$O_{ij}=\delta _{ij}[\underset{l=1}{\overset{\mathrm{}}{}}\gamma _l\overline{T}_{il}(w)T_{jl}(w)+\mathrm{c}.\mathrm{c}.].$$
(65)
### A Analytic Results

If the colloidal size $`a`$ is small, we can obtain an analytic expression. One may truncate the infinite series in eq.(65) at $`l=1`$:

$$O_{ij}\delta _{ij}\gamma _1\left[\overline{T}_{i1}(w)T_{j1}(w)+T_{i1}(w)\overline{T}_{j1}(w)\right].$$
(66)
However we notice that the $`(n+1)`$th and $`(n+2)`$th row of matrix $`𝐎^{(2)}`$ in eq.(64) can be written as

$`O_{(n+1)j}^{(2)}`$ $`=`$ $`\psi _j(w)={\displaystyle \frac{T_{j1}(w)}{\sqrt{\pi }}}`$ (67)
$`O_{(n+2)j}^{(2)}`$ $`=`$ $`\overline{\psi }_j(w)={\displaystyle \frac{\overline{T}_{j1}(w)}{\sqrt{\pi }}}.`$ (68)
If we multiply (67) by $`\sqrt{\pi }\gamma _1\overline{T}_{i1}(w)`$, multiply (68) by $`\sqrt{\pi }\gamma _1T_{i1}(w)`$, and add them together to the $`i`$th row, the off-diagonal part of $`O_{ij}`$ in eq.(66) is cancelled competely. Thus the final answer does not contain $`\gamma _1`$. This is true even if we do not make any truncations in eq.(66). Physically, it originates from the fact that the first angular momentum channel is occupied by the fixed ion, so it is not accessible to the free ions. Another way of saying this is that the fixed ion serves as a ‘Laughlin quasi hole’ in the quantum electron system.

To obtain a non-trivial result, we consider the second angular momentum channel $`l=2`$:

$$O_{ij}\delta _{ij}\gamma _2\left[\overline{T}_{i2}(w)T_{j2}(w)+T_{i2}(w)\overline{T}_{j2}(w)\right].$$
(69)
Applying Schur’s theorem, we have:

$`Z^{(2)}`$ $`=`$ $`det𝐎det\left(\left[\begin{array}{c}\mathrm{\Psi }(w)\\ \mathrm{\Psi }(\overline{w})\end{array}\right]𝐎^1\left[\begin{array}{cc}\mathrm{\Psi }^{}(w)& \mathrm{\Psi }^{}(\overline{w})\end{array}\right]\right)_{2\times 2}`$ (73)
$``$ $`Z_1Z_2.`$ (74)
We first concentrate on the second factor $`Z_2`$. It describes the two fixed ions (represented by wavefunctions $`\mathrm{\Psi }`$ and $`\mathrm{\Psi }^{}`$) interacting through the plasma medium (represented by matrix $`𝐎^1`$). For small $`\gamma _2`$, the matrix $`𝐎^1`$ can be approximated by:

$$O_{ij}^1\delta _{ij}+\gamma _2\left[\overline{T}_{i2}(w)T_{j2}(w)+T_{i2}(w)\overline{T}_{j2}(w)\right],$$
(75)
and $`Z_2`$ can be calculated to be (suppressing the $`(w)`$ arguments):

$`Z_2`$ $``$ $`\left[\left({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}|^2\right)^2\left|{\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}^2\right|^2\right]`$ (80)
$`+2\gamma _2\left[({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}|^2)\right|{\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}T_{k2}|^2`$
$`+\left({\displaystyle \underset{k=1}{\overset{n}{}}}|T_{k1}|^2\right)\left|{\displaystyle \underset{k=1}{\overset{n}{}}}\overline{T}_{k1}T_{k2}\right|^2`$
$`\left({\displaystyle \underset{k=1}{\overset{n}{}}}\overline{T}_{k1}^2\right)\left({\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}\overline{T}_{k2}\right)\left({\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}T_{k2}\right)`$
$`\left({\displaystyle \underset{k=1}{\overset{n}{}}}T_{k1}^2\right)\left({\displaystyle \underset{k=1}{\overset{n}{}}}\overline{T}_{k1}T_{k2}\right)\left({\displaystyle \underset{k=1}{\overset{n}{}}}\overline{T}_{k1}\overline{T}_{k2}\right)]`$
We can interprete the above result by thinking of the colloidal voids as a perturbation to the plasma system plus fixed ions. The first square bracket in eq.(80) is the result for the unperturbed system, namely, the ion-ion interaction via the plasma medium. It is precisely the two-body density function we obtained in eq.(31). The second square bracket in eq.(80) is the first order correction, namely, a fixed ion on one site interacts with the colloidal void on the other site. This ion-void interaction is usually found to be attractive, and since it involves higher order TC sums, its range is longer than that of the ion-ion repulsion (see results below). So far, we have missed another first order effect, namely, the fixed ion interacts with the colloidal void on the same site. The excess background charge underneath the void should reduce the effective ionic charge somehow. It turns out that this renormalization effect is captured precisely by the first factor $`Z_1`$ in eq.(74):

$$Z_1=det𝐎12\underset{k=1}{\overset{n}{}}|T_{k2}|^2\gamma _2.$$
(81)
In the bulk of plasma, we assume $`w=𝐢\frac{R}{2}`$ and use TC sum-rules (A20) and (A21) to obtain:

$$Z^{(2)}(12\gamma _2)\left(1\mathrm{e}^{R^2}\right)+2R^2\mathrm{e}^{R^2}\gamma _2.$$
(82)
In the above equation, the first term is the renormalized ion-ion repulsion and the second term is the ion-void attraction. The attraction is proportional to $`R^2\mathrm{e}^{R^2}`$, which has a longer range than the repulsion. This feature is not found in mean field approximations, thus must be caused by the strong correlation effects. Notice that the void-void repulsion and the renormalization of the ion-void attraction are both of second order, so they do not appear in eq.(82). The repulsion and attraction respectively dominates the short distance and long distance behavior, and the minimum free energy $`\beta F=\mathrm{ln}Z^{(2)}`$ occurs where the first derivative vanishes:

$$\frac{\mathrm{d}Z^{(2)}(R^{})}{\mathrm{d}R}=0\mathrm{o}rR^{}=\sqrt{\frac{1}{2\gamma _2}}.$$
(83)
Fig. 6 plots the calculated free energy vs. $`R`$ for various values of $`a`$. Note that when $`\gamma _2`$ is small ($`a=0.6`$), eq.(83) predicts that the location of the free energy minimum is far away from origin, in which case the attraction is strongly suppressed by the decaying factor $`\mathrm{e}^{R^2}`$, thus becoming virtually invisible.

Near the plasma edge, we use TC sum-rules (A30) through (A34) to obtain (for asymptotically large $`R`$):

$`Z^{(2)s}`$ $``$ $`\left[(A_{11}^s)^2{\displaystyle \frac{\mathrm{e}^{h^2}}{2\pi R^2}}\right]`$ (86)
$`+2\gamma _2[{\displaystyle \frac{(A_{11}^sh^2+A_{22}^s)\mathrm{e}^{h^2}\frac{2h}{\sqrt{2\pi }}\mathrm{e}^{\frac{3h^2}{2}}}{2\pi R^2}}`$
$`+A_{11}^s({\displaystyle \frac{\mathrm{e}^{h^2}}{\sqrt{2\pi }}}A_{11}^sA_{22}^s)],`$
Similarly for double edged system, we use TC sum-rules (A37) through (A41) to obtain

$`Z^{(2)d}`$ $``$ $`\left[(A_{11}^d)^2{\displaystyle \frac{\mathrm{e}^{h^2}}{2\pi R^2}}\mathrm{sin}^2(hR)\right]`$ (89)
$`+2\gamma _2\{{\displaystyle \frac{[A_{11}^dh^2\mathrm{cos}^2(hR)+A_{22}^d\mathrm{sin}^2(hR)]\mathrm{e}^{h^2}}{2\pi R^2}}`$
$`(A_{11}^d)^2A_{22}^d\},`$
### B Numerical Results

Our numerical results for charged colloids in the bulk case are plotted in Fig. 7. When $`a`$ is small ($`a=0.60.8`$) the numerical results agree very well with the analytic ones presented in Fig. 6. However, when $`a`$ is large ($`a=1.0`$), the two differ significantly in that the numerical free energy has an additional secondary maximum. This peak reflects the even longer ranged, second order effects which are not accounted for in the analytic studies.

We also tried putting double ion charges at the center of each colloids. This case is very similar to the divalent ions which are studied in Section III C. Our numerical results are presented in Fig. 8, where the colloidal sizes are chosen to be $`\sqrt{2}`$ times as large as those in Fig. 7, so that the two figures may be comparable. We see that the qualitative behavior is very similar. The attraction minimum stays further out.

Near the plasma edge, we use squares to represent colloidal particles as we did in Section IV B. The numerical results for two different sizes of colloids are presented in Fig. 9 and Fig. 10. In case I ($`a=0.8`$, Fig. 9), when the colloids are far away from the edge ($`h=5.0`$) they repel; as they become closer to the edge ($`h=2.3`$) they start to attract each other within a certain range ($`R<2.25`$), beyond that range they still repel; at a particular point ($`h=1.5`$), the range of the attraction seems to extend to infinity; finally, when they are too close to the edge ($`h=0.8`$), due to insufficient screening, a long range repulsion is found. In case II ($`a=0.9`$, Fig. 10), the colloids are larger in size and attract each other even in the bulk ($`h=5.0`$); the attraction becomes deeper ($`h=2.0`$) and longer ranged ($`h=1.3`$) when the colloids move towards the edge.

Fig. 11 shows the free energy of colloids in a double edged system. When the two boundaries are widely separated ($`h=5.0`$), the colloids repel; as the separation decreases, the free energy develops a long oscillating tail. These are consistent with the analytic results. Note that the oscillation ‘wavelength’, according to eq.(89), is set by the strip width $`h`$, and is independent of the size $`a`$ of the colloids. This conclusion is confirmed by Fig. 12.

Finally, Fig. 13 shows the results for doubly charged colloids near a single plasma edge. They are qualitatively similar to those shown in Fig. 10.

### C Where Does Colloid Stay?

To better understand the charged colloidal interactions, especially near the plasma edge, it is instructive to learn where they want to stay in the plasma. We put a charged colloid at $`w=r\mathrm{e}^{𝐢\theta }`$, which is $`\frac{h}{2}`$ away from the plasma edge at $`(r+\frac{h}{2})`$, and calculate the partition function to be:

$`Z^{(1)s}`$ $`=`$ $`det𝐎^{(1)}`$ (90)
$`𝐎^{(1)}`$ $`=`$ $`\left[\begin{array}{cc}𝐎& \mathrm{\Psi }^{}(w)\\ \mathrm{\Psi }(w)& 0\end{array}\right],`$ (93)
where the $`n\times n`$ overlap matrix $`𝐎`$ is:

$$O_{ij}=\delta _{ij}\underset{l=1}{\overset{\mathrm{}}{}}\gamma _l\overline{T}_{il}(w)T_{jl}(w).$$
(94)
It can be shown that the first angular momentum channel is blocked by the fixed ion charge so that the final answer does not depend on $`\gamma _1`$. For small $`a`$, we truncate the infinite summation in eq.(94) at $`l=2`$ and use Schur’s theorem to find:

$`Z^{(1)s}`$ $``$ $`\left[1{\displaystyle \underset{i=1}{\overset{n}{}}}|T_{i2}|^2\gamma _2\right]\left[{\displaystyle \underset{i=1}{\overset{n}{}}}|T_{i1}|^2+\left|{\displaystyle \underset{i=1}{\overset{n}{}}}\overline{T}_{i1}T_{i2}\right|^2\gamma _2\right]`$ (95)
$`=`$ $`A_{11}^s+\left({\displaystyle \frac{\mathrm{e}^{h^2}}{2\pi }}A_{11}^sA_{22}^s\right)\gamma _2.`$ (96)
Fig. 14 shows the numerical results done for a square charged colloid. We plot the free energy $`\beta F`$ as a function of $`h`$ for three different sized squares. In all cases, we find the colloid is attracted to a free energy minimum near the plasma edge (at $`h=0`$). For $`h<0`$, the free energy quickly blows up. This tells that the colloid does not want to flee away from the plasma edge, because its fixed charge is attracted by the background charge and wants to be in the lowest available angular momentum channel. The optimal distance $`h^{}`$ is in the range in which attraction between colloids occurs.

## VI Summary

In this paper, we considered a strongly coupled ($`\beta =2`$) 2D one-component plasma model and calculated the effective interaction between two colloids both in the bulk and near boundaries. We found that due to ion-void interaction, the colloids may attract despite the fact that ion-ion and void-void interactions are both repulsive. Near the plasma edge, the attraction effects become more prominent. Though our model is oversimplified, the exact results provide useful insights into the problem, and can serve as a testing ground for more sophisticated approximate treatments.

## Acknowledgement

We thank Adrian Parsegian and Rudi Podgornik for introducing us to this problem. RR is grateful to SMG and the Physics Department of Indiana University for supporting a summer visit during which this work was initiated. SMG and RR are grateful to the Institute for Theoretical Physics at UCSB where part of this work was performed. This work is supported by NSF DMR-9714055.

## A Properties of Lowest Landau Level Wavefunctions

Some important properties of lowest Landau level (LLL) wavefunctions are summarized here. They are frequently used throughout the paper.

### 1 Orthonormality and Completeness

The LLL wavefunctions defined in eq.(5) are orthonormal:

$$i|j=dz\overline{\psi }_i(z)\psi _j(z)=\delta _{ij}.$$
(A1)
Sometimes, we want to calculate a partial inner product (the overlap integrated over part of the entire plane), and the result is generally not diagonal. However, if the integration region has azimuthal symmetry, the conservation of angular momentum guarantees an exception. As an example, let us take the region to be a circle with radius $`a`$ centered at origin (denoted by the subscript),

$$i|j_0=_0^a\rho d\rho _0^{2\pi }d\varphi \frac{\rho ^{i+j2}\mathrm{e}^{\rho ^2}\mathrm{e}^{𝐢(ji)\varphi }}{\pi \sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(j)}}=\delta _{ij}\gamma _i,$$
(A2)
where the imcomplete gamma function $`\gamma _i`$ is defined as the following:

$$\gamma _i\frac{1}{\mathrm{\Gamma }(i)}_0^{a^2}du\mathrm{e}^uu^{i1}=1\mathrm{e}^{a^2}\left(\underset{k=0}{\overset{i1}{}}\frac{a^{2k}}{k!}\right).$$
(A3)
Notice that in the limit $`a0`$, $`\gamma _iO(a^{2i})`$.

The completeness tells that any state lying in LLL can be expressed as a linear combination of these basis wavefunctions. An example is the wavefunction $`\psi _i(z)`$ with its argument translated by a vector $`w`$:

$$\mathrm{e}^{\frac{1}{2}(\overline{z}wz\overline{w})}\psi _i(z+w)=\underset{j=1}{\overset{\mathrm{}}{}}T_{ij}(w)\psi _j(z).$$
(A4)
Here, the extra pure phase factor in front of $`\psi _i(z+w)`$ is needed to project the translated wavefunction into the LLL.

### 2 Translation Coefficients

The linear coefficients $`T_{ij}`$s defined in eq.(A4) are termed translation coefficients (TCs). Setting $`z=0`$ in eq.(A4), we find:

$$\psi _i(w)=\underset{j=1}{\overset{\mathrm{}}{}}T_{ij}(w)\psi _j(0)=\frac{T_{i1}(w)}{\sqrt{\pi }}.$$
(A5)
Thus the LLL basis functions are just special cases of TCs. Another useful application of TC is to calculate the partial inner product on a circle with radius $`a`$ centered at $`w0`$ (denoted by the subscript):

$`i|j_w`$ $`=`$ $`{\displaystyle _{|z|<a}}dz\overline{\psi }_i(z+w)\psi _j(z+w)`$ (A6)
$`=`$ $`{\displaystyle _{|z|<a}}dz{\displaystyle \underset{k,l}{}}\overline{T}_{ik}(w)\overline{\psi }_k(z)\psi _l(z)T_{jl}(w)`$ (A7)
$`=`$ $`{\displaystyle \underset{k,l}{}}\overline{T}_{ik}(w)k|l_0T_{jl}(w)`$ (A8)
$`=`$ $`{\displaystyle \underset{l=1}{\overset{\mathrm{}}{}}}\gamma _l\overline{T}_{il}(w)T_{jl}(w).`$ (A9)
The TCs can be evaluated explicitly. We multiply $`\overline{\psi }_k(z)`$ to both sides of eq.(A4), and integrate $`z`$ over the entire plane. On the right hand side, this picks out one particular $`T_{ik}`$ from the summation because of the orthonormality:

$`T_{ik}(w)`$ $`=`$ $`{\displaystyle dz\mathrm{e}^{\frac{1}{2}(\overline{z}wz\overline{w})}\overline{\psi }_k(z)\psi _i(z+w)}`$ (A10)
$`=`$ $`{\displaystyle \frac{\mathrm{e}^{|w|^2/2}}{\pi \sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(k)}}}{\displaystyle dz\mathrm{e}^{|z|^2z\overline{w}}\overline{z}^{k1}(z+w)^{i1}}`$ (A11)
$`=`$ $`{\displaystyle \frac{\mathrm{e}^{|w|^2/2}}{\pi \sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(k)}}}{\displaystyle dz\mathrm{e}^{|z|^2}}`$ (A13)
$`\times \left[{\displaystyle \underset{t=0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{(\overline{w})^tz^t}{t!}}\right]\overline{z}^{k1}\left[{\displaystyle \underset{s=1}{\overset{i}{}}}C_{s1}^{i1}z^{s1}w^{is}\right]`$
$`=`$ $`\sqrt{\mathrm{\Gamma }(i)\mathrm{\Gamma }(k)}(1)^k\widehat{w}^{ik}\mathrm{e}^{\frac{|w|^2}{2}}`$ (A15)
$`\times {\displaystyle \underset{s=1}{\overset{\mathrm{min}(i,k)}{}}}{\displaystyle \frac{(1)^s|w|^{i+k2s}}{\mathrm{\Gamma }(s)(is)!(ks)!}}.`$
Where, in the 3rd equation we used Taylor and binomial expansions, and in the 4th equation we used orthonormality (A1). Eq.(A15) immediately tells that the TCs have the following inversion symmetry:

$$T_{ij}(w)=\overline{T}_{ji}(w)=T_{ji}(\overline{w}).$$
(A16)
Another important property of TC arises when one composes two successive translations:

$$T_{ij}(w_1+w_2)=\mathrm{e}^{\frac{1}{2}(\overline{w}_1w_2w_1\overline{w}_2)}\underset{k=1}{\overset{\mathrm{}}{}}T_{ik}(w_1)T_{kj}(w_2).$$
(A17)
Or using eq.(A16), this becomes:

$`T_{ij}(w_1+w_2)`$ $`=`$ $`\mathrm{e}^{\frac{1}{2}(\overline{w}_1w_2w_1\overline{w}_2)}{\displaystyle \underset{k=1}{\overset{\mathrm{}}{}}}\overline{T}_{ki}(w_1)T_{kj}(w_2)`$ (A18)
$`=`$ $`\mathrm{e}^{\frac{1}{2}(\overline{w}_1w_2w_1\overline{w}_2)}{\displaystyle \underset{k=1}{\overset{\mathrm{}}{}}}T_{ki}(\overline{w}_1)T_{kj}(w_2)`$ (A19)
As a simple application, we set $`w_1=w_2=w`$ in eq.(A18):

$$\underset{n\mathrm{}}{lim}\underset{k=1}{\overset{n}{}}\overline{T}_{ki}(w)T_{kj}(w)=T_{ij}(0)=\delta _{ij},$$
(A20)
or set $`w_1=w_2=\frac{𝐢R}{2}`$ in eq.(A19):

$$\underset{n\mathrm{}}{lim}\underset{k=1}{\overset{n}{}}T_{ki}\left(\frac{𝐢R}{2}\right)T_{kj}\left(\frac{𝐢R}{2}\right)=T_{ij}(𝐢R).$$
(A21)
Note that if we interprete $`T_{ij}`$ as the matrix elements of $`𝐓`$, then eq.(A20) has a simple matrix form:

$$𝐓^{}(w)𝐓(w)=1,$$
(A22)
showing that matrix $`𝐓`$ is unitary. In the same matrix notation, eq.(A16) is rewritten as

$$𝐓^{}(w)=𝐓(w)$$
(A23)
which along with eq.(A22) makes sense since $`𝐓(w)`$ should be the inverse of $`𝐓(w)`$.

Eq.(A20) and (A21) are examples of TC sum-rules. They are frequently used in analysing the bulk of the plasma. Near the edge, we often encounter TC summations of a different type:

$`A_{ij}^s`$ $``$ $`\underset{r\mathrm{}}{lim}{\displaystyle \underset{k=1}{\overset{(r+\frac{h}{2})^2}{}}}\overline{T}_{ki}(r\mathrm{e}^{𝐢\theta })T_{kj}(r\mathrm{e}^{𝐢\theta })`$ (A24)
$`B_{ij}^s`$ $``$ $`\underset{r\mathrm{}}{lim}{\displaystyle \underset{k=1}{\overset{(r+\frac{h}{2})^2}{}}}T_{ki}(r\mathrm{e}^{𝐢\theta })T_{kj}(r\mathrm{e}^{𝐢\theta }).`$ (A25)
Here, the superscript ‘s’ denotes single edged system and $`\theta \frac{R}{2r}`$. The difference between the bulk TC sum and edge TC sum is as following: treated as a perturbation, colloids in the bulk mainly perturb the low angular momentum channels. Such a perturbation has a discrete characteristic, namely only a small number of channels are strongly affected. Colloids near the boundary, on the other hand, perturb the high angular momentum channels. The perturbation has a continuum nature, namely there are a large number of channels ($`r^2hr<k<r^2+hr`$) being weakly affected. As a result, the evaluation of edge TC sum inevitably requires one to convert the discrete summation into a continuum integral. For example, let us calculate $`A_{11}^s`$:

$`A_{11}^s`$ $`=`$ $`\underset{r\mathrm{}}{lim}{\displaystyle \underset{k=1}{\overset{(r+\frac{h}{2})^2}{}}}{\displaystyle \frac{r^{2k2}}{\mathrm{\Gamma }(k)}}\mathrm{e}^{r^2}`$ (A26)
$``$ $`\underset{r\mathrm{}}{lim}{\displaystyle _0^{r^2+hr}}{\displaystyle \frac{\mathrm{d}x}{\sqrt{2\pi x}}}\mathrm{e}^{2x\mathrm{ln}r+xx\mathrm{ln}xr^2},`$ (A27)
where we have defined $`xk1`$ and used Stirling’s asymptotic formula:

$$\mathrm{\Gamma }(x+1)\sqrt{2\pi x}\mathrm{e}^{x+x\mathrm{ln}x}.$$
(A28)
The integal in eq.(A27) can be evaluated using the saddle point approximation, in which one expands the exponent near its minimum at $`x=r^2`$ up to the second order:

$$A_{11}^s\underset{r\mathrm{}}{lim}_0^{r^2+hr}\frac{\mathrm{d}x}{\sqrt{2\pi }r}\mathrm{e}^{\frac{1}{2r^2}(xr^2)^2},$$
(A29)
substituting $`y\frac{xr^2}{r}`$, one finds

$$A_{11}^s=_{\mathrm{}}^h\frac{\mathrm{d}y}{\sqrt{2\pi }}\mathrm{e}^{y^2/2}=\frac{1}{2}\left[1+\mathrm{\Phi }\left(\frac{h}{\sqrt{2}}\right)\right].$$
(A30)
Where $`\mathrm{\Phi }(x)`$ is the standard error function. Other TC sum rules can be obtained through quite similar procedures. We omit the tedious derivations and simply list the results below:

$`A_{22}^s`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[1+\mathrm{\Phi }\left({\displaystyle \frac{h}{\sqrt{2}}}\right){\displaystyle \frac{h\mathrm{e}^{h^2/2}}{\sqrt{2\pi }}}\right],`$ (A31)
$`A_{12}^s`$ $`=`$ $`{\displaystyle \frac{\mathrm{e}^{h^2/2}}{\sqrt{2\pi }}},`$ (A32)
$`B_{11}^s`$ $`=`$ $`{\displaystyle \frac{\mathrm{e}^{h^2/2}}{\sqrt{2\pi }R}}\mathrm{e}^{𝐢(hR\frac{\pi }{2})},`$ (A33)
$`B_{12}^s`$ $`=`$ $`{\displaystyle \frac{h\mathrm{e}^{h^2/2}}{\sqrt{2\pi }R}}\mathrm{e}^{𝐢(hR\frac{\pi }{2})}.`$ (A34)
Note that the above results are correct only asymptotically for $`R\mathrm{}`$ ($`R2r\theta `$).

In a double edged system, we encounter the same type of TC summation, except that the lower limit is changed from $`k=1`$ to $`k=(r\frac{h}{2})^2`$:

$`A_{ij}^d`$ $``$ $`\underset{r\mathrm{}}{lim}{\displaystyle \underset{k=(r\frac{h}{2})^2}{\overset{(r+\frac{h}{2})^2}{}}}\overline{T}_{ki}(r\mathrm{e}^{𝐢\theta })T_{kj}(r\mathrm{e}^{𝐢\theta })`$ (A35)
$`B_{ij}^d`$ $``$ $`\underset{r\mathrm{}}{lim}{\displaystyle \underset{k=(r\frac{h}{2})^2}{\overset{(r+\frac{h}{2})^2}{}}}T_{ki}(r\mathrm{e}^{𝐢\theta })T_{kj}(r\mathrm{e}^{𝐢\theta }).`$ (A36)
Correspondingly, we have: (the superscript ‘d’ denotes double edges):

$`A_{11}^d`$ $`=`$ $`\mathrm{\Phi }\left({\displaystyle \frac{h}{\sqrt{2}}}\right),`$ (A37)
$`A_{22}^d`$ $`=`$ $`\mathrm{\Phi }\left({\displaystyle \frac{h}{\sqrt{2}}}\right){\displaystyle \frac{h\mathrm{e}^{h^2/2}}{\sqrt{2\pi }}},`$ (A38)
$`A_{12}^d`$ $`=`$ $`0,`$ (A39)
$`B_{11}^d`$ $`=`$ $`{\displaystyle \frac{\mathrm{e}^{h^2/2}}{\sqrt{2\pi }R}}\mathrm{sin}(hR)\mathrm{e}^{𝐢(hR\frac{\pi }{2})},`$ (A40)
$`B_{12}^d`$ $`=`$ $`{\displaystyle \frac{h\mathrm{e}^{h^2/2}}{\sqrt{2\pi }R}}\mathrm{cos}(hR)\mathrm{e}^{𝐢(hR\frac{\pi }{2})}.`$ (A41)
## B Schur’s Theorem

For any four matrices $`𝐀_{r\times r}`$, $`𝐁_{r\times s}`$, $`𝐂_{s\times r}`$ and $`𝐃_{s\times s}`$, if $`𝐀`$ is not singular, Schur’s theorem states

$$det\left[\begin{array}{cc}𝐀& 𝐁\\ 𝐂& 𝐃\end{array}\right]=det𝐀det(𝐃\mathrm{𝐂𝐀}^1𝐁).$$
(B1)
The proof is straightforward. We consider the following matrix identity:

$$\left[\begin{array}{cc}1_{r\times r}& 0_{r\times s}\\ \mathrm{𝐂𝐀}^1& 1_{s\times s}\end{array}\right]\left[\begin{array}{cc}𝐀& 𝐁\\ 𝐂& 𝐃\end{array}\right]=\left[\begin{array}{cc}𝐀& 𝐁\\ 0_{s\times r}& 𝐃\mathrm{𝐂𝐀}^1𝐁\end{array}\right].$$
(B2)
Taking the determinant of both sides of eq.(B2), eq.(B1) follows immediately.

Schur’s theorem has a very useful corollary: for any matrices $`𝐉_{r\times s}`$ and $`𝐊_{s\times r}`$, one has

$$det(1\mathrm{𝐉𝐊})_{r\times r}=det(1\mathrm{𝐊𝐉})_{s\times s}.$$
(B3)
We can prove this via very similar method. First we consider

$$\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 0_{s\times r}& 1_{s\times s}\end{array}\right]\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 𝐊& 1_{s\times s}\end{array}\right]=\left[\begin{array}{cc}(1\mathrm{𝐉𝐊})_{r\times r}& 0_{r\times s}\\ 𝐊& 1_{s\times s}\end{array}\right],$$
(B4)
which tells

$$det\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 𝐊& 1_{s\times s}\end{array}\right]=det(1\mathrm{𝐉𝐊})_{r\times r};$$
(B5)
on the other hand, if we consider

$$\left[\begin{array}{cc}1_{r\times r}& 0_{r\times s}\\ 𝐊& 1_{s\times s}\end{array}\right]\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 𝐊& 1_{s\times s}\end{array}\right]=\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 0_{s\times r}& 1\mathrm{𝐊𝐉}\end{array}\right],$$
(B6)
we obtain

$$det\left[\begin{array}{cc}1_{r\times r}& 𝐉\\ 𝐊& 1_{s\times s}\end{array}\right]=det(1\mathrm{𝐊𝐉})_{s\times s}.$$
(B7)
A combination of eq.(B5) and (B7) results in eq.(B3). Letting $`𝐉=\gamma _1𝐊^{}`$ in eq.(B3), we obtain eq.(52).