GROMACS version:
GROMACS modification: Yes/No
Here post your question
Hi, Sir/Ma’am
I want to include three types of phosphoinositide lipids (PIP1, PIP2, PIP3) in my simulation. However, during the minimization step, I encountered an error (attached screenshot) indicating that the molecule type POP1 is not found (I have given even phosphoinositide.itp file separately).
Can someone help me to address this error?
I would suggest not modifying the itp file. Change your topology file instead. But make sure that you are using the correct molecule type. I can’t know if you should use POP1_3, POP1_4 or POP1_5.
I followed your suggestions and modified the topology to perform the minimization, but I still encountered the “total potential energy is infinite” error along with an atom mismatch warning (attached screenshot). I ran the simulation for 40 fs and 20 fs, and I also adjusted the emtol value to 0.001. However, I’m consistently getting LINCS errors during the minimization step.
I reran my simulation from scratch but encountered the same error again. To neutralize my system, I used the genion command, but it didn’t resolve the issue. When I run the genion command, it replaces one of my lipid molecules DPG3.
What steps can I take to resolve this error? I have attached my command line and the corresponding file for your reference. systemtop.txt (1.1 KB) minimizationmdp.txt (4.6 KB)
There’s no need trying to run the simulation until you get gmx genion to neutralize it properly. How are you selecting the solvent group? Are you using an index file?
Going through attached images in semi-random order is not very easy. Could you please copy and paste all input and output from the command-line as text instead? Preferably to an attached text file.
Never get used to use the -maxwarn option in gmx grompp unless you really know what you are doing. It might help you get simulations started, but the results might not be correct. -maxwarn 8 is a sign that there are lots of problems. Always go through each step of the preparation carefully and verify that the system and setup is correct.
You also seem to use options to neutralize the system already with insane-new-pip.py. Isn’t the system neutral after that? If not, why not? Why would you need gmx genion at all?
I increased the salt concentration from 0.1 to 0.5 in an insane step, which resulted in the “system has non-zero total charge” error again during the minimization step. After running the genion command, the error was resolved, but I encountered an infinite total potential energy error during minimization, along with the replacement of the DPG3 lipid molecule in the neutralized system.gro file. In the genion step, I chose to replace water molecules (W) with NaCl ions. I also set -maxwarn 8 as a precaution, although I did not face any issues other than the aforementioned error and an atom mismatch warning. I’ll ensure not to use -maxwarn unless absolutely necessary. Also, I’ll send the text file you mentioned soon.
It’s difficult to say. From the insane output it says that the system has a charge of -154. Is that before or after adding ions? I presume it’s before, is that correct? Because in the gmx make_ndx output there is 154 more Na+ than Cl-, so presumably the system should be neutral at that stage.
Then you run gmx grompp to generate ions.tpr. At that point you have a net charge of -557 it should be 0. You shouldn’t even have to run gmx genion. That’s the point in your protocol when you need to take a step back. Have a look at your topology file. Why is there a net charge? Are the itp files for the components different than the charge that insane assumes? Does any component have the wrong charge? Is the net charge of 7wc4 the same in insane and in the itp file? The topology file you posted before does not match the console output from gmx make_ndx, e.g.:
24 NA+ : 7435 atoms
25 CL- : 7281 atoms
compared to:
NA+ 7360
CL- 7206
The net charge of -154 is neutralized in both cases, though.
After reviewing your comments, I re-prepared my protein using CHARMM-GUI and completed a 10 ns equilibration. I then took that structure for the MARTINI simulation, and fortunately, I did not encounter any issues with atom overlap or infinite potential energy errors. However, I did face a different problem during the minimization step, specifically with atom 485. I would appreciate any suggestions you might have to solve this issue; I’ve attached the error screenshot and text file for your reference.
I attempted to resolve this by changing the time step (dt) from 0.001 to 0.02, adjusting the energy threshold (emtol) values to 0.001, 1, and 10, and setting constraints to none, but the minimization did not complete successfully.
Additionally, I have been having problems with genion, as it removes my DPG3 lipid molecule from the .gro file. My system was not neutralized, and I am still receiving warnings about having a non-zero total charge. As you mentioned, the system were neutral up to the insane stage, but it appears to become more negatively charged during the grompp step. I’m currently cross-referencing both my insane script and lipid .itp file to identify any discrepancies.
I’ll keep you updated once I find the issue.
Thanks a lot for your help and suggestions
It is not necessarily a problem that the energy minimization does not converge. If the following NVT equilibration works OK, then I would not worry too much. If the equilibration crashes, then you would have to look at your coordinates to see why there is such a large force on atom 485.
As I said before, if the system is already neutral (and even with a specifically set ion concentration) after generating it, there is no reason to consider using gmx genion at all. If gmx grompp complains about the system not being neutral, you will have to go through your itp files one at a time, comparing the charges to their expected charges, and see where the net charge comes from.
If you can’t find the source of the problem, simplify your setup. Start by using only one lipid type in the membrane. When that is working (you can run gmx grompp without any warnings at all, or notices about the system being net charged), you can add one more lipid type. That way, you should be able to find when the problem occurs.
I am trying to add Gb3 lipid to my membrane system. Since this lipid is not present in the insane script, I modified the DPG3 lipid by editing its ITP file, removing extra atoms to create a new Gb3 lipid ITP. I renamed DPG3 to Gb3 in the GRO file accordingly.
However, during the minimization step, I encountered a coordinate mismatch error, which indicates that my system.top file needs to be updated to include 1314 additional atoms. I used a command (taken from your old post
) to calculate the number of lipid molecules, and here are the current counts in my system.gro file:
POPC (Lipid): 8773
POPE (Lipid): 5257
POPI (Lipid): 2731
DPSM (Lipid): 2014
CHOL (Lipid): 14149
Gb3 (Lipid): 3797
SAP5 (Lipid): 307
SAP6 (Lipid): 325
SAP3 (Lipid): 685
W (Water): 348390
NA+ (Sodium ion): 5231
CL- (Chloride ion): 5053
I also attached my system.top file, but I am struggling to correctly update the atom types and lipid counts in the system.top file to fix the mismatch error. Could you help guide me on how to resolve this issue?
The topology file you posted has a completely different count. Are you presenting the number of atoms in each molecule type in the list above? I’m sorry, but I won’t go through all your molecule types to find out how many atoms there are in each. I don’t mean to be unhelpful, but I’m sure you can find ways to count the number of molecules you’ve got of each type in your gro file and update the topology file yourself. If you can’t think of anything else then:
Start with the first molecule type in your coordinate file (7wc4) count the number of atoms of that molecule type (i.e., scroll until you find the next molecule type).
Check the number of atoms in the corresponding itp file (or in the topology file itself if it’s written there).
Divide the number you counted in 1 by the number you counted in 2. For 7wc4 I expect you get 1 molecule.
Make sure that the number you got in 3 matches the number of the corresponding molecule type (in the order they are listed) in the topology file.
Redo from 1 with the next molecule type in your coordinate file.
Keep in mind that all molecules of the same type do not necessarily follow each other. E.g., have a look at POPC in your topology file. I’m not saying this is the quickest way to do it, but it will always work, or you are doing something wrong.
You should make sure that you have an itp that matches exactly what you want to simulate - it must match the number of atoms in the coordinate file. Modifying an itp file by just removing some atoms is risky, unless you really know what you are doing and verify that the partial charges are still correct. If there is not a matching itp file you will have to consider how to handle that - and be ready to explain it to a potential reviewer.