DPD of water beads


I am new to doing DPD. I am trying to run DPD simulation of water beads in LAMMPS using real units.
But density of my simulation box isn’t coming close to 1. But for water if my parameters are correct, shouldn’t it come close to 1.

Each bead consists 6 water molecules, m = 108.09 g/mol, and degree of coarse-graining (Nm)= 6
I have calculated DPD repulsion parameter using, a = Kb.T.(K.Nm-1)/(2.alpha.rho_dpd) …from [1]
rho_dpd = 3
k = 15.9836 (water compressibility at 300K)
aplha = 0.101
I am getting, a = Kb.T(156.60)
At 300K and on per mol basis I am getting a = 93.351 kcal/mol
However, documentation states that “a” units should be in force, i.e., (kcal/mol)/Angstroms. So, my question is with what length am I supposed to divide my “a” value to get it in these units?
I tried dividing by rcut, which in my case is 8.14 Angstroms, but when I ran NPT simulation at 300K and 1atm, my density keeps on decreasing.

Also, documentation states that gamma units has to be force/velocity. But in [2], Steve suggests to take gamma unit as fs.(kcal/mol), which is dimensionally not equal to force/velocity. Instead, shouldn’t it be fs(kcal/mol)/Angstroms^2?
I am using gamma as 0.09 (g/mol)/fs.

May be it is something else that I missed while doing DPD simulations, if someone could provide some insight I would be really grateful.

Lastly, I am providing my LAMMPS script, and the thermo output as well.

I don’t know DPD well enough to comment. But there is a large literature on running DPD

models with LAMMPS. I suggest you look for papers similar to your model and contact
the authors if you are still perplexed.