I want to perform an energy minimization of a complex protein/ligand/crystal water. I created the topology files, included the ligand, etc. Everything seemed to work fine. However, when I compared the position of some water molecules after the minimization with the crystal, they changed a lot. Most of the oxygen atoms changed their position slightly, but some appeared on the other side of the protein. Is there a way I could hinder the water molecules from changing so much (be optimized, but not appear on the other side)?
I read in the manual that I could turn PBC off and optimize the system in a “vacuum” by adding the PBC = no, nstlist, and cuttofs = 0, but there’s the message “But vacuum simulations are (temporarily) not supported”. So I don’t know if it is possible.
If you read this thread: How to treat specific water molecules as ligand? you will get some advice how to do it. In brief, you need to make a new molecule type for restrained water molecules and you need to replace settles with constraints in that molecule type. You also need to add a restraints section (see other restrained molecules). Then in the [ molecules ] section in your topology you need to make sure that the crystal waters are represented by the restrained water molecule type (you have to count how many water molecules come before), e.g.:
[ molecules ]
protein 1
SOL 1001
SOL_restrained 1
SOL 234
SOL_restrained 1
SOL 123456
Then you have to consider if you need to restrain them only during minimization or also during equilibration. If you need to restrain them during production as well, chances are that there is a high rate of exchange (quite probable) or that the crystal water positions are not very favourable in the simulations.