Genion with alchemical transformations

User discussions
1

Hello GROMACS Community!

I’m currently engaged in alchemical transformation simulations involving changes in net charge. I’m familiar with the challenges associated with charge changes and have been using the dual-system-in-a-single-box method to address them. However, during some test runs, I’ve come across a perplexing issue: In my test system, the net charge differs between state A and state B. For initial simulations, I run each state individually, either at lambda=0 for state A or lambda=1 for state B.

Here’s what puzzles me: In both lambda states, I use gmx genion to neutralize the system. Strangely, the number of ions added by genion remains consistent across both lambda states, aligning only with the neutralization needs of state A. This implies that for state B simulations, the number of counter ions isn’t appropriate.

So, my question is: How exactly does genion determine the net charge for neutralization? Does it consider the lambda state specified in the mdp file to accurately gauge the system’s net charge, or does it default to using the lambda=0 state for this calculation?

Your insights on this would be greatly appreciated.

Best regards,
FullmetalAlchemist

2

Hi, if you’re using alchemical ions to ensure charge neutrality when (un-)charging alchemical residues, you should create another [ moleculetype ] (using a normal ion as a template) that will alchemically change the charge. AFAIK gmx genion doesn’t do that at all.

If you have, say, 50 K molecules and you create a K_alch alchemical molecule, you can edit your [ molecules ] section to write:

[ molecules ]
...
K 49
K_alch 1

K_alch can also have the same atom type and name as in K, only the name in [ moleculetype ] and the B-state charge should be different.

For a general discussion of the problem, also check section 7.2.2. in this excellent review.