XX (and more generally XY) spin models became very popular lately, especially in connection to quantum-information theoretic problems. The XY model (without magnetic field though) has been introduced and properly solved in 1961 by <PERSON>, <PERSON> and <PERSON> (Ann. Phys. 16, 407, (1961)). By the way this is a great reference to start with, if you didn't read this paper already. Nonetheless, in many recent articles, this model is solved incorrectly, and this might be the origin of the discrepancy that you observed. Let me address your point.
First I imagine that the average that you mean is the ground state average: i.e. you digonalize the Hamiltonian, take the eigenvector related to the lowest eigenvalue, and compute statistical averages with this vector. This said, to specify completely the problem you must assign some boundary conditions (BC). For many reasons (which I won't enter here), the most appropriate boundary condition is the periodic one (PBC) for the spins. With this I mean $\sigma^\alpha_{L+1}=\sigma^\alpha_1$, $L$ being the system size. 976
When you map the spin problem to a Fermionic one through the Jordan-Wigner (JW) mapping, the Hamiltonian becomes that of free Fermions plus a peculiar border term. This is because the JW string accumulates in the bond which couples spins 1 and L. The result is that an XY (or XX) spin model with PBC is equivalent to a free-fermion model with a parity dependent boundary conditions (the parity operator is $P\sim \prod \sigma^z_i \sim (-1)^N$, and $N$ is the total number operator in the Fermi picture, and the $\sim$ symbol ignores potential factors of $i$). In the parity 1 sector (which comprises the ground state) BC are antiperiodic (i.e. $c_{L+1}=-c_1$ for the Fermi operators $c_i$). This boils down to antiperiodic quasimomenta when you do the Fourier transform: $k_n = \pi (2n+1)/L$ with $n=0,1,\ldots,L-1$. If you instead choose open boundary conditions, there is no funky border term in the fermionic problem either, but the diagonalization of an open chain is slightly more cumbersome.
This could be one source of errors.
Another potential source of error is degeneracy (but this shouldn't be your case). When you diagonalize the Hamiltonian you should check that the ground state is unique. In a quasi-free systems you have degenerate ground states when the one-particle dispersion becomes exactly zero for some quasimomentum. This for instance, is known to happen in a pure free fermionic chain with nearest neighbour interaction and PBC, when the size $L$ is a multiple of 4.
As I said, a great reference is <PERSON>, <PERSON> and <PERSON> article, where the precise diagonalization procedure is presented. Boundary conditions are also discussed in a recent paper of mine (PRA 81, 060101 (2010)) or in the book of <PERSON>, Conformal invariance and critical phenomena.
Hope this helps.
Best, <PERSON> 976