\sqrt{B_0}+ \hat O^{x}_2+ \hat W \}$, by grouping two operators into one $\hat W=\hat O^{x}_3+\hat O^{x}_4$. Then, from bosonization, it follows that $\sum_r S^{y}_{r}S^{y}_{r+1} \to  \int \mathrm{d} r  \{ \sqrt{B_0} - \hat O^{x}_2+ \hat W \}$. For $\lambda \to 1$, using Eq. (\[ee1\]) for the case of $a=x$ and $b=y$, we obtain $3\langle\langle \hat W(0) \hat W(r)\rangle \rangle =\langle\langle \hat O^{x}_2(0) \hat O^{x}_2(r)\rangle \rangle$ and thus, $$\label{mainresult}
G^{a,a}_u(r)=\frac{4}{3}G^{x,x}_{B_2}(r)=\frac{1}{6\pi^4}  \frac{\ln^2r}{ r^4}.$$ This exact asymptotic expression is our main result.

![Bilinear spin correlation function of the Heisenberg spin-$\frac{1}{2}$ AFM chain $(-1)^r\rho G^{a,a}_r\cdot 10^3$, where we have introduced the cord distance on the circle $\rho= L/\pi\sin(r\pi/L)$. Data shown are for $L=96$ spins with periodic boundary conditions (PBC) and are obtained by keeping typically $m=1500$ states in DMRG simulations. Bullets indicate numerical data [@comment1] and continuous line is the curve, $10^3d/\ln^{\frac{3}{2}}(c\rho)$, with $d\simeq 0.067$ and $c\simeq 16$ coefficients obtained by fits of analytical curve to numerical data for $r>10$. Inset shows $(-1)^r\rho G^{a,a}_r\cdot 10^3$, for $L=144$ spins and for distances, $40 \le r\le 72$. Fitting to the data for $L=144$ sites gives similar estimates for $d$ and $c$.[]{data-label="fig:pd"}](Fig1.eps){width="7.5cm"}

Numerical results
=================

In the remaining part we will present a numerical check of Eq. (\[mainresult\]) 