GROMACS version: 2020.3

GROMACS modification: No

I am performing Brownian Dynamics with hard spheres. But the simulation stops after a random number of steps, Gromacs keeps running but I don’t have evolution on the number of steps in any output file.

First I was using only the C12 part of the LJ potential and I thought this problem could be due to two particles getting to close together and the potential getting to big. So I added the C6 term, the problem took longer to appear, the simulation was able to run for a longer time, but it is still there.

I am putting below the mdp file I am using, the itp for the particles and the topology.

The system is 1000 particles in a box of 5 x 5 x 5 um3.

Any help is greatly appreciated! Thank you so much!

*mdp file*_

ld-seed = 55560 ; Use random seed

;################################# INTEGRATOR ##################################

integrator = bd ; Langevin thermostat

dt = 1000 ; Timestep (ps)

nsteps = 20000000 ; Simulation duration (timesteps)

nstcomm = 0 ; Center of mass motion removal interval

comm-mode = None ; Center of mass motion removal mode

bd_fric = 24000 ; The friction coefficient for each particle is calculated as mass/tau-t.

;################################## ENSEMBLE ###################################

ref_t = 295 ; System temperature (K)

tau_t = 2.0 ; Thermostat time constant (ps)

tc_grps = system ; Apply thermostat to complete system

;############################## IMPLICIT SOLVENT ###############################

implicit_solvent = no ; NO implicit solvent

;########################### NONBONDED INTERACTIONS ############################

cutoff_scheme = Verlet ; Method of managing neighbor lists

pbc = xyz ; Periodic Boundary Conditions in all directions

coulombtype = cut-off ; Calculate coulomb interactions using PME

rcoulomb = 10 ; Coulomb cutoff

vdwtype = cut-off ; Calculate van der Waals interactions using PME

rvdw = 10 ; Van der Waals cutoff

rlist = 10 ; Neighbor list cutoff

nstlist = 100 ;

;################################### OUTPUT ####################################

nstlog = 5000 ; Log output interval (timesteps)

nstenergy = 5000 ; Energy output interval (timesteps)

nstcalcenergy = 5000 ; Energy calculation interval (timesteps)

nstxout = 1000 ; Trajectory output interval (timesteps)

nstvout = 5000 ; Velocity outout interval (timesteps)

nstfout = 5000 ; Force output interval (timesteps)

nstxout-compressed = 0 ;

*itp file*_________

; Topology of solid particles.

; No include forcefield parameters

[ defaults ]

; nb-type combination-rule gen-pairs fudgeLJ fudgeQQ

1 1 no

[ atomtypes ]

;type mass charge ptype c6 c12

FA 4.8e+04 0.000 A 3.898e+04 6.090e+08

FB 3.0e+08 0.000 A 2.494e+12 2.494e+24

FC 1.03e+3 0.000 A 1.818e+03 1.326e+06

[nonbond_params]

FA FA 1 3.898e+04 6.090e+08

FB FB 1 2.494e+12 2.495e+24

FC FC 1 1.818e+03 1.326e+06

[ moleculetype ]

;name nrexcl

F005 1 ; sigma=05nm Epsilon= 0.6236 kJmol-1 mass= 48000 Da C12= 6.08984e+08 kJmol-1nm12

[ atoms ]

; nr type resnr residu atom cgnr charge mass

1 FA 1 F005 FA 0 0.00 4.8e+04

[ moleculetype ]

;name nrexcl

F100 1 ; sigma=100nm Epsilon= 0.6236 kJmol-1 mass= 3.0e+08 Da C12= 2.49440e+24 kJmol-1nm12

[ atoms ]

; nr type resnr residu atom cgnr charge mass

1 FB 1 F100 FB 0 0.00 3.0e+08

[ moleculetype ]

;name nrexcl

F003 1 ; sigma=03nm Epsilon= 0.6236 kJmol-1 mass= 1029 Da C12= 1.32563e+06 kJmol-1nm12

[ atoms ]

; nr type resnr residu atom cgnr charge mass

1 FC 1 F003 FC 0 0.00 1.02914e+03

*top file*___

; Topology of 1000 particles.

#include “solidParticlesC6.itp”

[ system ]

; systemname

1000 particles 005nm

[ molecules ]

; molecule number

F005 1000