5
$\begingroup$

I need to compute the helicity modulus as a function of temperature for a three-dimensional X-Y model (see N.K. Kultanov, Yu.E. Lozovik, "The critical behavior of the 3D X-Y model and its relation with fractal properties of the vortex excitations"). I am resorting to using a Monte Carlo simulation because I cannot find a closed expression for the superfluid density in a 3-D MC model in the literature. I have written a simple Monte Carlo simulation and I am experiencing some problems, i.e. long thermalization times and noisy results even when averaging on a huge number of simulations.

Since my code is simple, are there any tricks I could use to make the computation more efficient? References and source code would be very much appreciated. I am primarily interested in the 3-D case, since my code seems to work well for the 2-D case.

If anyone is interested I've uploaded my code (in C) here: http://pastebin.com/ScrMEaQ8

Improve this question
$\endgroup$
9
  • 2
    $\begingroup$ zakk, Monte Carlo is a programming technique. Making a code run faster is a programming questions. There are some related physics questions about what you can ignore or elide for your purposes but you haven't asked those. You asked "How can I make this program faster?" And you haven't even told us what you've tried along those lines. $\endgroup$
    – dmckee --- ex-moderator kitten
    Jun 30, 2012 at 17:02
  • 2
    $\begingroup$ @dmckee I have to disagree. Monte Carlo is not a programming technique, it's a general framework for exploring statistical distributions. Granted people don't do MC calculations by hand (any more) but OP isn't just asking us to vectorize his loops, there is a canonical answer to this question. $\endgroup$
    – wsc
    Jun 30, 2012 at 17:31
  • 1
    $\begingroup$ I know I'm on shaky ground invoking this logic, but both of the papers I cite below were seen as sufficiently important general, conceptual leaps to justify publishing in PRL (as opposed to say PRE, or Comp. Phys. Comm.) and are both papers that all budding young mechanics-of-the-statistical-persuasion should know -- even if they never write any codes! $\endgroup$
    – wsc
    Jun 30, 2012 at 17:52
  • 1
    $\begingroup$ @zakk: Could you include in your question a description of the problem you are solving, as well as the methods and algorithms you used to solve it? It would help us out tremendously to give you better feedback on how to make it better. $\endgroup$
    – Paul
    Jun 30, 2012 at 18:02
  • 1
    $\begingroup$ @zakk: Huge statistical errors are often a function of being near a critical point. You'll need to provide some more solid information about the specifics of your XY model in order to get more specific answers, just like Paul suggested. $\endgroup$
    – aeismail
    Jul 1, 2012 at 15:35

2 Answers 2

Reset to default
9
$\begingroup$

The long thermalization time that you're running into is a generic problem that typically goes under the name "critical slowing down" and is common to the local-update scheme that you're using (you update by locally changing a single spin at a time).

Once you realize that, the way out is to do better sampling - local updates are out so you have to invent global updates. Two great ways of doing this are as follows:

1) Cluster updates using the Wolff algorithm:

http://prl.aps.org/abstract/PRL/v62/i4/p361_1

(Sorry, no free version available)

2) Find a "dual" model to the one you're studying, with the same partition function but new degrees of freedom. For Ising and XY models if you play this game you end up with "bond-current" models where the currents are updated globally using a worm algorithm detailed nicely in

http://arxiv.org/abs/cond-mat/0103146v1

Oh! And as I googled for that Prokofiev reference, I see that someone even implemented it with cute graphics in Mathematica!

http://demonstrations.wolfram.com/WormAlgorithmForJCurrentModel/

Improve this answer
$\endgroup$
0
3
$\begingroup$

Just taking a quick look at your code, there are a few things that could be improved.

  • You use three random number calls to pick a site, when you could just use modulo arithmetic with one number from 1 to $N = N_x \times N_y \times N_z$, and compute the values of $i$, $j$, and $k$.
  • Similarly, in update_configuration, you are performing loops that lead to a function call; you'd probably be better off if you reorganized the function to contain the loops instead. That way you can take better advantage of memory structures and data parallelism than you have now.
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.