8
$\begingroup$

I have written some molecular dynamics code that I want to use to model known chemical reactions.

By known I mean that it has been observed in a lab that $\ce{A +B->C}$.

I am not interested in finding out if it is feasible that two things react (i.e. is there a favorable energy configuration?), rather I am interested in finding out what the concentration over time looks like given a volume, concentration, temperature and a pair of chemicals that I know will react.

At the moment I just check if they are very close to each other. This produces reasonable looking concentration curves.

  • Is this is a correct way of doing things? I am looking for publications so it has to be defend-able, not just OK.

I am also interested in modelling large bio-molecules in a coarse grained simulation.

  • Are there any particular potential functions that I should use for large biomolecules, or will LJ do?
Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ Usually MD only simulations can't handle reactions, it is some hybrid of QM and MD. Alternatively you could consider a reaction probability given an intermolecular distance (i.e., if the molecules are inside a reaction radius, grab a random uniform number. If it is > .6 then reagent molecules are removed and a product is put in its place.) You'd have to be careful with artifacts in your velocity and acceleration terms in the MD simulation with a heavy handed technique like that. MD/QM is probably the best bet but computationally expensive. $\endgroup$
    – scs217
    Commented Aug 9, 2013 at 19:02
  • $\begingroup$ What MD package are you using by the way? $\endgroup$
    – scs217
    Commented Aug 9, 2013 at 19:02
  • $\begingroup$ I wrote my own package, as i'm trying to simulate a very specific thing. If I matched the reaction probabilities to published results do you think that would work? I am going to be simulating well understood reactions $\endgroup$
    – RNs_Ghost
    Commented Aug 9, 2013 at 19:10
  • $\begingroup$ I suppose it depends on what you're trying to get out of the MD simulation. If you're just looking for a concentration curve, MD is overkill. If you want to have a complete time history and see where the reactions are taking place, then MD is a good option. I worry about conservation of momentum if you're just popping molecules in and out and how that might introduce artifacts into the simulation. $\endgroup$
    – scs217
    Commented Aug 22, 2013 at 14:12
  • 2
    $\begingroup$ If you want to publish this thing, you may want to be a little more familiar with existing solutions, literature etc. Your question suggest you try to re-invent the wheel which can be fun hobby, but a waste of time if it is your research. $\endgroup$
    – Greg
    Commented Sep 27, 2015 at 22:38

1 Answer 1

Reset to default
2
$\begingroup$

Since you've written your own code, you will probably be able to judge this best, but I'll give you some pointers as to where I'd start.

As far as I've gathered, you simulate a NVT ensemble of a pair of chemicals.

First of all you have to decide whether you want explicit representation of the solvent or not. Hint: The latter is computationally far more expensive and I doubt that you have implemented all the optimizations that you'd find in larger MD software.

This is because this choice can drastically influence the diffusion properties of your chemicals. In a coarse-grained force field they might still be able to diffuse around pretty much unhindered, while when using explicit atomistic representation of the solvent, the properties may change to a larger degree.

The way I'd do it would be to check whether the two areas of the molecules that react with each other are close to one another, meaning that some distance criteria $(d<d_\text{rxn.})$ has to apply. When this criteria is met, you can make the reaction happen by creating new bonds and changing the atom types that are bonded on-the-fly. It wouldn't be so tricky to implement this part, I guess. Caveat: This will only work if you have a reaction whose mechanism is close to the proposed scheme. You might introduce nontrivial errors otherwise.

Then you can model the "real" concentration curves of all the substances, and simultaneously change your system to reflect the change in number of molecules. I wouldn't be surprised that you'd find that the pressure would slowly decrease as the reaction takes place (since you keep the volume constant).

Improve this answer
$\endgroup$
2
  • $\begingroup$ Would it not be even better, if the reaction parameters are known and the transition state is known, to put a criterion on a function of distance and velocity along the reaction coordinate? This might allow to include also the temperature parameter within the MD simulation. The other problem with the distance approach is that it would not allow that there could also be a reverse reaction (equilibrium situation), whereas one could also check the distance / velocity criterion along the dissociation coordinate of C. $\endgroup$
    – PLD
    Commented May 15, 2014 at 14:12
  • $\begingroup$ @PLD If the thermalization algorithm is good enough, you would already have enlarged diffusion and thus a higher rate of contact between A and B, so the temperature is already factored in here. Additionally, you could simulate the reverse reaction by assigning a probability to the formation of the bond. Or even assign a probability to the breaking of the bond, so to the reverse reaction directly. $\endgroup$
    – tschoppi
    Commented May 15, 2014 at 15:26

Your Answer

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

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