Fix rigid/npt: Non-numeric box dimension - simulation unstable

Hello,

I am attempting to simulate a box of CO2 with fix rigid/nvt. When I simulate it in in nvt, the pressure of the system is very negative, so I have been trying to first use fix rigid/npt to manage the pressure of the system.

When it runs the fix rigid/npt though I get the error “Non-numeric box dimension - simulation unstable.” I have tried this on multiple different systems, and even when the pressure is already at a normal, reasonable value and I attempt to fix it at the same pressure that it was at before I get this error. I am using the create_atoms command to randomly put the linear CO2 molecules in a simple box, so I also am not sure how this could be caused by non-physical inputs.

I have pasted my input file and the molecule file I am using to create the CO2 molecules below. I would appreciate any help with this situation.

Thanks,
Michael

input file:

TraPPE Carbon Dioxide NVT

units real
atom_style full
boundary p p p
region box block 0 46.39721813432450 0 46.3972181343245 0 46.3972181343245
create_box 2 box bond/types 1 angle/types 0 &
extra/bond/per/atom 2 extra/angle/per/atom 1 extra/special/per/atom 2

mass 1 12.0107 #C
mass 2 15.9994 #O

variable dt equal dt
variable set_temp equal 290.00
variable Tdamp equal 100.0

#290K, 10 MPa, 19.951 mol/L
pair_style lj/cut/coul/long 8.5
pair_modify tail yes
bond_style harmonic
angle_style harmonic
kspace_style pppm 1.0e-4

pair_coeff 1 1 0.053654828567 2.8000
pair_coeff 2 2 0.156970053954 3.05000

dihedral_style none
improper_style none

bond_coeff 1 5000.00 1.160

molecule co2 TraPPE.mol
create_atoms 0 random 1200 213487 NULL mol co2 25367 overlap 1.37

group carbon type 1
group oxygen type 2
set group carbon charge 0.7
set group oxygen charge -0.35

velocity all create 290.00 432467 rot yes dist gaussian

reset_timestep 0

fix 1 all rigid/npt molecule temp {set_temp} {set_temp} ${Tdamp} iso 98.6923 98.6923 10.0
thermo 1000
run 500000
unfix 1
reset_timestep 0

molecule file:

Carbon Dioxide molecule. TraPPE geometry

3 atoms
2 bonds
0 angles

Coords

1 0.00000 0.00000 0.00000
2 1.16300 0.00000 0.00000
3 -1.16300 0.00000 0.00000

Types

1 1 # C
2 2 # O
3 2 # O

Charges

1 0.700
2 -0.350
3 -0.350

Bonds

1 1 1 2
2 1 1 3

1 Like

I notice two issues in your input file:

  1. 10 (fs in real unit) is too small as the relaxation time of the barostat (as your timestep is, by default of the real unit, 1 fs). After changing the iso 98.6923 98.6923 10.0 to iso 98.6923 98.6923 1000 it does not throw that error anymore. You probably copied the examples on the doc page, though they’re likely more reasonable for other units settings.
  2. in create_atoms 0 random 1200 213487 NULL mol co2 25367 overlap 1.37, the 1.37 (≈the minimum distance between atoms of different CO2 molecules in the initial structure) is likely too small. If you include Fmax in thermo output you’ll find that the maximum force applied on an atom is ~40000 real units on the first step, which is shockingly large (common values are on the order of 10^1 or 10^2). That can easily lead to crashing simulations/unphysical results, as the molecules are gaining too much velocities/temperature in a single step.
3 Likes

I have made these changes and the simulation works a lot better now, thanks for your help.