GROMACS version:
GROMACS modification: Yes/No
Here post your question
Hi! I am trying to equilibrate my membrane using the C-rescale barostat with three different rounds of relaxation decreasing position restraints from 1000 to 500 to 100 before production run. However, I’m noticing that there’s very little difference between the membranes between each relaxation period and it looks very crystalline. In addition, when running the production run with the Parrinello-Rahman barostat, it crashes due to LINCS warnings.
I’m thinking that it might be because the tau-p = 2 paired with Parrinello-Rahman’s instability with fluctuations might’ve resulted in the production crashing so I might increase that or just try C-rescale for production run as that’s worked for me in the past but I’m not sure why the membrane is hardly moving and relaxing in between equillibrations. Any advice or suggestions would be lovely. Thank you in advance!!
The issue is that my membrane looks very crystalline in between equillibrations (NPT 1ns each, NVT is 100ps). I’m thinking this is might be contributing to the crashes in the production run, but I’m not totally sure. I guess I’m expecting the membrane to relax more and display more movement compared to a very symmetrical and rigid starting state for the production run.
in that case you need to increase the time for equilibration. Though, I would suggest run atleast a micro second simulation (production )and check in between like 100 ns or 200 ns etc.
You are starting from an extremely well organised bilayer (I guess you geometrically built it by equi-spacing some lipid structure file?)
It is very tricky to equilibrate a system like that. You membrane is still looking very crystal-like in the images because 100 kJ/mol/nm² is still quite a large restraining potential. If you are applying it to all the atoms of the lipids then yes, it makes sense that they are barely relaxing. Moreover, I feel like you do not have enough water molecules, as the water slabs on the two sides seem very very thin, but maybe once the box relaxes on the xy plane they will become thicker.
You should for sure run longer relaxation runs with weaker potentials. Also, you may start applying restraints only to specific parts of the systems, i.e. a restraint on a couple of atoms between the lipid heads and the lipid tails only along the z axis. This should i) let your lipids move laterally, changing the packing of the system and letting the box relax also on xy plane ii) let the heads reorient a bit and relax and iii) let the tails disrupt the well packed structure and go towards a more physiological gel/liquid state (depending on the lipid types and the temperature).
You can also consider to run high temperature to help the later diffusion of lipids and try to break down as much as possible the initial crystal correlation and packing. Eventually, you will have to check the thickness, area per lipid, order parameters etc and compare with some exp result/some comp paper where they use your same force field, just to be sure that your system has properly relaxed and is at equilibrium.
Just a question: can’t you use CHARMM-GUI to generate the starting structure?