8
$\begingroup$

I have a polymer and solvent on a lattice, and I am performing an NVT Monte Carlo simulation. I am performing the following moves on it:

enter image description here

My question is, how should these moves be conducted to effectively sample phase space? I have 4 functions written up: end_rotation(), kink_jump, crankshaft(), reptation(), which I can call sequentially in my code - but I want to understand what is the best way to call them to sample as much of this ensemble as I can to get meaningful statistics from this - something like the Flory exponent.

$\endgroup$
3
  • $\begingroup$ Do you know of any references on this topic (for example, I assume the image in the post came from such paper)? I just want to avoid a case where a user points you toward a potential reference only to find out you had already read papers X, Y, and Z. $\endgroup$
    – Tyberius
    Nov 7, 2021 at 1:23
  • $\begingroup$ Also, I'm not an expert in this area, but I would suspect a decent starting point would just be to randomly select a type of move at each step. If some type of move is clearly being tried more often than needed, you can adjust the distribution away from being a uniform selection. $\endgroup$
    – Tyberius
    Nov 7, 2021 at 2:09
  • $\begingroup$ Hello @Tyberius. This image came from one of Scott Shell's lectures at UCSB sites.engineering.ucsb.edu/~shell/che210d/… $\endgroup$
    – megamence
    Nov 14, 2021 at 23:03

1 Answer 1

7
$\begingroup$

This is a very good and tricky question, which I don't think has a clear and definite answer. I think I should also preface by saying that I can't answer it from the point of view of polymer physics, but instead from my experience in biomolecular simulations. However, this particular problem is basically a version of the sampling exploration-exploitation trade-off and as such is equivalent in many different fields.

You can think about two types of sampling problems: one with a small discrete number of states where your aim is to explore all of them during the course of your simulation, or alternatively (and much more commonly) a setting where you have too many states to explore and you just want to explore as many as possible.

In the former case, it is conceptually (although not practically) straightforward to define a clear-cut "best" way to sample based on transition matrices: in this case, you can assign an effective averaged transition matrix (with dimensions equal to the number of discrete states) to each of the different moves. One such example is a simulated tempering simulation with a small number of discrete temperatures, where your aim is to have maximum mobility between the lowest and the highest temperature. This defines an optimisable quantity which is a function of these transition matrices and you can in this way determine the best type of move (if all moves are ergodic), or the best combination of moves otherwise. It is not obvious to me how you can actually do the latter in practice, but I imagine it is possible to do adaptively using some sort of reinforcement learning technique, since you have a countable action space and a well-defined reward function. Definitely not trivial, but well-defined and possible. This would make for a nice PhD project!

The second setting (i.e. the one in question) is much harder, because the objective is not even well-defined. On one hand you want to explore as many states as possible, but on the other hand you want to explore all of them "evenly". What I mean by this is that even if one of the moves was guaranteed to provide you with strictly new configurations for the duration of your whole simulation, you would still want to attempt the other moves, since they provide configurations that cannot be explored by the first one (e.g. if you only did successful end rotations, the sampling would still be pretty useless). Therefore, you need to ask yourself how much you care about each of these moves and assign each one a "target" weight, which is purely arbitrary and based on your goals. The simplest and most obvious way to do this is to give an equal importance to each of these moves.

But that's not the end of it, because in practice each of these moves will have different average acceptance probabilities, making some of them more efficient than the others. Again, if you were simply interested in new configurations regardless of the type of move, you should be simply able to adaptively determine the one giving you the best acceptance and choose it randomly proportional to its acceptance (how you define this proportionality is also of course arbitrary, because the only way to be completely sure certain moves are bad is to evaluate them and you don't want to do that!). However, if you care about all the moves to an extent, then you need to define how much performance loss you are willing to tolerate for each one of them. This is also completely arbitrary and it is not even clear to me if you can have a objectively good way of defining that.

So what should you do? One way is to simply attempt all the moves with equal probability and accept that this is not the most efficient way, but is also the easiest way to define. Another (probably better) way is to continuously estimate the "usefulness" of each move and attempt it with a probability directly proportional to its currently estimated average acceptance: e.g. if one move has a 100 times better average acceptance than another one, it will be attempted 100 times as frequently. I have seen a variant of this idea be utilised in the case of constant-pH simulations, where you want to change as many protonation states as possible, but you only want to do that if you are reasonably sure that the protonation state is likely at the pH of interest.

TL;DR: It is completely arbitrary how you do that, but it might make sense to start with equal probabilities and at some point switch to probabilities that are proportional to the average estimated acceptance of each of the moves.

$\endgroup$
6
  • 2
    $\begingroup$ +10. Welcome back Godzilla!!! Master of molecular dynamics!!! I hope to see you around even more often, and thanks for all your contributions!!! $\endgroup$ Nov 14, 2021 at 4:09
  • 2
    $\begingroup$ Thanks Nike! I have been quite busy lately but I try to monitor the MD questions every once in a while and answer if I can be useful :) $\endgroup$
    – Godzilla
    Nov 14, 2021 at 7:59
  • 1
    $\begingroup$ This is a very informative answer. Thank you! I am trying to write a detailed balance, but getting the numbers down for that balance is, as you mentioned, not at all clearly defined. $\endgroup$
    – megamence
    Nov 14, 2021 at 23:05
  • $\begingroup$ @megamence I completely ignored the topic of balance, but yes, any of these moves need to satisfy detailed balance (or just balance more generally). If you are doing adaptation (as I suggested), you are by definition breaking Markovianity, but as long as the quantities you are adaptively modifying are asymptotically convergent, it is still a valid sampling algorithm. Thus, if you base your proposal probabilities on the success of every move so far, this quantity will eventually converge to the "true" efficiency, and while non-Markovian, this is a valid sampling algorithm. $\endgroup$
    – Godzilla
    Nov 14, 2021 at 23:45
  • 2
    $\begingroup$ I would add that the relative rates (expected or observed) of the moves could factor into the sampling especially if there is a goal to approximate something physically meaningful. So if the slowest (least probable) move is sampled a sufficient number of times, then all the others are as well. $\endgroup$ Nov 15, 2021 at 5:46

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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