Im doing a Metropolis Monte Carlo simulation with particles on a sphere and have a question concerning the random movement in a given time step.
I understand that to obtain a uniform distribution of random points on the sphere to begin with it is not enough to use the naive simplest way (use spherical coordinates with constant R and pick random angles theta and phi), but one has to for example use one of these methods: http://mathworld.wolfram.com/SpherePointPicking.html
Looking at another code for a Monte Carlo on a sphere I see a fairly complicated way to generate random moves: pick a random particle, calculate the rotation matrix moving it to the north pole, find a random cartesian vector less than a certain length and move it to the north pole, normalize the cartesian vector and then rotate it back to the vicinity of the original particle position.
This is all to get an unbiased new random position. I don`t understand the rationale completely although I suspect it is connected to detailed balance. But my feeling is that there should be an easier (i.e. faster) way to do this. Actually, intuitively I feel that in this case it is ok to find two small random angles theta and phi and add them to the position of the particle – or would this be a mistake?
Best Answer
$\newcommand{\+}{^{\dagger}}% \newcommand{\angles}[1]{\left\langle #1 \right\rangle}% \newcommand{\braces}[1]{\left\lbrace #1 \right\rbrace}% \newcommand{\bracks}[1]{\left\lbrack #1 \right\rbrack}% \newcommand{\ceil}[1]{\,\left\lceil #1 \right\rceil\,}% \newcommand{\dd}{{\rm d}}% \newcommand{\ds}[1]{\displaystyle{#1}}% \newcommand{\equalby}[1]{{#1 \atop {= \atop \vphantom{\huge A}}}}% \newcommand{\expo}[1]{\,{\rm e}^{#1}\,}% \newcommand{\fermi}{\,{\rm f}}% \newcommand{\floor}[1]{\,\left\lfloor #1 \right\rfloor\,}% \newcommand{\half}{{1 \over 2}}% \newcommand{\ic}{{\rm i}}% \newcommand{\iff}{\Longleftrightarrow} \newcommand{\imp}{\Longrightarrow}% \newcommand{\isdiv}{\,\left.\right\vert\,}% \newcommand{\ket}[1]{\left\vert #1\right\rangle}% \newcommand{\ol}[1]{\overline{#1}}% \newcommand{\pars}[1]{\left( #1 \right)}% \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\pp}{{\cal P}}% \newcommand{\root}[2][]{\,\sqrt[#1]{\,#2\,}\,}% \newcommand{\sech}{\,{\rm sech}}% \newcommand{\sgn}{\,{\rm sgn}}% \newcommand{\totald}[3][]{\frac{{\rm d}^{#1} #2}{{\rm d} #3^{#1}}} \newcommand{\ul}[1]{\underline{#1}}% \newcommand{\verts}[1]{\left\vert\, #1 \,\right\vert}$ Let's assume we have a mechanism to generate random numbers in $\left[0,1\right)$ and ${\rm P}\pars{\Omega_{\vec{r}}}$ is a distribution function for random points in a sphere of radius $a > 0$. $\Omega_{\vec{r}}$ is the solid angle. In this case, ${\rm P}\pars{\Omega_{\vec{r}}}$ is, indeed, $\Omega_{\vec{r}}\,\,$-independent: $$ 1 = \int{\rm P}\pars{\Omega_{\vec{r}}}\,\dd\Omega_{\vec{r}} = {\rm P}\pars{\Omega_{\vec{r}}}\int\dd\Omega_{\vec{r}} = {\rm P}\pars{\vec{r}}\pars{4\pi} \quad\imp\quad{\rm P}\pars{\vec{r}} = {1 \over 4\pi} $$ Then, $$ 1 = \int_{0}^{\pi}\half\,\sin\pars{\theta}\,\dd\theta\int_{0}^{2\pi} \,{\dd\phi \over 2\pi} $$ We can generate random numbers $\xi_{\theta}$ and $\xi_{\phi}$ such that: $$ \bracks{~\half\,\sin\pars{\theta}\,\dd\theta = \dd\xi_{\theta}\,, \quad\xi_{0} = 0 \imp \theta = 0~}\ \mbox{and}\ \bracks{~{\dd\phi \over 2\pi} = \dd\xi_{\phi}\,,\quad\xi_{0} = 0 \imp \phi = 0~} $$ Those relations yield: $\ds{\half\bracks{-\cos\pars{\theta} + 1} = \xi_{\theta}}$ $\ds{\pars{~\mbox{or/and}\ \sin\pars{\theta/2} = \root{\xi_{\theta}}~}}$ and $\ds{\phi = 2\pi\,\xi_{\phi}}$: $$\color{#0000ff}{\large% \theta = 2\arcsin\pars{\root{\xi_{\theta}}}\,,\qquad \phi = 2\pi\xi_{\phi}} $$ As we mentioned above, $\xi_{\theta}$ and $\xi_{\phi}$ are uniformly distributed in $\left[0, 1\right)$.