Dear GROMACS users, I have a question for calculating the potential of mean force between 2 plates (i.g., graphene) in water. It is not technically as simple as it seem. I found this question has been asked many times (I sum them up in the bottom of this note, so maybe it will be helpful for someone with similar question in future) but I cannot see a clear answer completely solving the problem. Information in the manual is quite limited on this topic as well.
There are two ways I can see to do the calculation: 1. umbrella sampling; 2. free energy perturbation (JACS, 2005, 127, 3556). Both of them, however, need to deal with 2 questions first in order to get any meaningful result: 1. keep the plates to be really plates (flat surface); 2. inhibit plate rotation. There are 3 ways that come into my mind to solve these 2 problems, but none of them works as I tried. (Similar idea has been suggested in the previous posts regarding to the problem.) I am using a simple NPT system in the free energy perturbation scheme, because it is physically more meaningful and simpler to test than umbrella potential. But problem comes from the 3 method I tried: 1) Merely constructing a really complicated itp file with bond, angle, dihedral, imp dihedral cannot easily maintain a planar plate. Rotation of the plates is not easy to deal with, which I still don't have any idea yet. (Someone suggested other methods rather than this one in the previous post.) 2) Freezing all plates will make the plates do what they should do, but there is problem with the pressure coupling according to the manual. GMX 4.0.x wouldn't work as specified in the manual. GMX 3.3.3 works as it seems but the system is very easy to blow up. 3) Use position RESTRAINT on all atoms on the plates in XYZ direction (basically fixing them in the initial position) as suggested in previous posts. As an example: ----------------------------------------------------------------------------------------------- [ position_restraints ] ; ai funct fc 1 1 95000 95000 95000 2 1 95000 95000 95000 3 1 95000 95000 95000 4 1 95000 95000 95000 .... every atom in the plates ------------------------------------------------------------------------------------------------ However, or small force constant, the plate won't be maintained. For large force constant (the above example), the system will easily blow up even with time step of 1fs. For sure, once we have too large vibration in the system, an even smaller time step is required. Then the efficiency is going to be really really low. (To provide further details: I minimize the system with both steep and cg first and run the 400 ps NVT(Nose) run before the NPT (Parrinello, isotropic or semiisotropic) run. The test system is fairly small: 570 waters (SPC/E) with 2 plates of graphene (60 atoms each, amber parameters). So there should not be problem with the equilibration. I tuned the force constant in the restraint to see whether the plate change or not.) I will really appreciate it if someone can suggest some way out or where my problem is in the simulation. Is it the right way of restraining the system? Is there better way of doing this? It will be great if discussions coming along. Thank you in advanced for your help! (There has been tons of paper doing this with some other package, I am sure GROMACS can do it. I am just not the person who can figure this out.) Thanks, Zhe Links for the previous representative post related: http://lists.gromacs.org/pipermail/gmx-users/2010-February/048615.html http://lists.gromacs.org/pipermail/gmx-users/2010-January/048326.html http://lists.gromacs.org/pipermail/gmx-users/2009-December/047696.html -- gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
