It works!, but I noticed that the Correlation (or Conditional Probabilities) $Pr(X{_i}=1|X{_k}=1)$ between any fixed pair (i and k) slightly change as we add more objects into the physical space. So, I changed $Cov(Y{_i},Y{_k})$ to be $w(||x{_i}-x{_k}||)$ and it worked with fixed conditional probabilities, even if I add more objects. I am still testing it, but does this make sense to you? For your information, my distance functions is $e^{(-\frac{eucaliandistance}{dcorr})}$

OK, I think I understand your suggestion now. I will test it right away.

Acually, It is positive definite, but does not have a dichotomized Gaussian distribution for the correlation matrix.

If i understood you correctly, what you are saying is to simply use the Cov(Yi,Yk) above as the covariance matrix and plug it in the bernoulli generator. If so, then Unfortunately, I tried that but I keep getting an Unacceptable correlation matrix, which seams to be not positive definite.