Hi! I am trying to use tabulated interactions in my molecular dynamics runs, but I am not sure I am doing it right.
I have a system of graphine sheet and methane molecules.
The molecule type for these two is GRA and MET, respectively. Since both of them consist of single atoms, the atom types of GRA and MET are C1 and CH4, respectively.
I have three tables in my working directory:
table.xvg (C1 C1 interactions)
table_C1_CH4.xvg
table_CH4_CH4.xvg
in my .mdp file I have included:
vdwtype = user
coulombtype = Cut-off
energygrps = C1 CH4
energygrp-table = CH4 CH4 C1 CH4
So, the question is whether my implementation makes sense.
Well, the issue is with my test particle insertion simulations (tpi), where the chemical potential gets -inf values at some point. I insert a CH4 molecule in a box of my graphene layer and 200 CH4 molecules. Please see below
Then I though to check just the CH4 molecules alone. I prepared a box of CH4 molecules and ran an NPT simulation to see whether I could reproduce the experimental density at given P and T, which worked fine. (308K and 20bar: experimental density = 12.8 kg/m3; my simulation density = 13.2 kg/m3)
Then I thought to try to do a tpi simulation in this bulk CH4 system (by inserting a CH4 molecule). I did and I again my chemical potential reached these -inf values:
Ah, but TPI with user tables is a bit of a special case as very short distances can occur. What have you done with the tables for distances close to zero? You need to ensure that the net potential, dispersion+repulsion+Coulomb, is always large and positive there.
First line I made by hand, to ensure at distance=0 no funny things go on. But maybe what you said doesnt reflect to my table’s small distances, and I need to make further adjustment?
That looks reasonable. The repulsion values at zero are quite close to max float. I wouldn’t think this can give issues, but it would be good to check if reducing that avoids the issue.
