Precise asymptotics are probably beyond the reach of current technology. Presently, all known methods to obtain lower bounds on the least eigenvalue (or singular value) of square random matrices require either an explicit formula for the joint density function (which is only available for very special ensembles, such as the Wishart ensemble), as done for instance in a paper of Edelman, or rely heavily on independence on the entries, as done in the recent papers of Rudelson-Vershynin or myself and Vu. Your ensemble does not seem to have enough algebraic structure for an explicit description of the density function, and also does not have enough independence to use the methods of Rudelson-Vershynin or Van Vu and myself.

But one can get upper bounds on the least eigenvalue by testing the matrix against some test vectors. For instance, one expects (from the birthday paradox) that two of the $x_i, x_{i'}$ should be at a distance $O(1/n^2)$ from each other. The two corresponding rows of the matrix then also differ by $O(1/n^2)$ when viewed one entry at a time, and so this already gives an upper bound of $O( n^{-3/2} )$ for the least eigenvalue. Perhaps one can do a bit better than this by locating multiple rows where the $x_i$ are close to each other and then using Taylor expansion to view the matrix as a perturbation of a low rank matrix, which could give better upper bounds. It would be difficult though to show that these bounds are sharp.

I think I can get an upper bound of $O(1/n^2)$ by exhibiting a vector $v$ of magnitude comparable to $1$ which gets mapped to a vector of magnitude $O(1/n^2)$. The basic idea is to exploit the birthday paradox to find (with high probability) two indices $i \neq i'$ such that $x_i-x_{i'} = O(1/n^2)$. It should also be possible to then find another additional index $i''$ such that $x_{i''} = x_i + O(1/n)$.

Now look at the $i^{th}$ and $(i')^{th}$ rows, which have components $1/(1+|x_i-x_j|)$ and $1/(1+|x_{i'}-x_j|)$. These rows differ by $O(n^{-2})$ in each coefficient. This already gives an upper bound of $O(n^{-3/2})$ for the smallest eigenvalue, but one can do better by using Taylor expansion to note that the difference between the two components is $(x_i-x_{i'}) \hbox{sgn}( x_i-x_j ) / (1 + |x_i-x_j|)^2 + O(n^{-4})$ except when $x_j$ is very close to $x_i$, at which point we only have the crude bound of $O(n^{-2})$. Similarly, the difference between the $i''$ and $i$ rows is something like $(x_i-x_{i''}) \hbox{sgn}( x_i-x_j ) / (1 + |x_i-x_j|)^2 + O(n^{-4})$ except when $x_j$ is too close to $x_i$. So we can use a multiple of the second difference to mostly cancel off the first difference, and end up with a linear combination of three rows in which most entries have size $O(n^{-4})$ and only about $O(1)$ entries have size $O(n^{-2})$. This seems to give an upper bound of $O(n^{-2})$ for the least eigenvalue (or least singular value), though I have not fully checked the details.

To get a matching lower bound is trickier. One may have to move to a Fourier representation of the matrix as this would more readily capture the positive definiteness of the matrix (as suggested by Bochner's theorem).