Hi all,

When I run LAMMPS using QEQ for ionic liquid systems with 200 ion-pairs, I always encounter this kind of error after several ns. Can you suggest how to deal with this problem? My input file is:

```
# Number of atoms = 9600; Total charge = 6.953437825529818e-12
units real
boundary p p p
atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style opls
special_bonds lj/coul 0.0 0.0 0.5
# remove hybrid if not necessary
pair_style lj/cut/coul/long 12.0 12.0
pair_modify mix geometric tail yes
kspace_style ewald 1.0e-5
read_data data-box.lmp
pair_coeff 1 1 0.170000 3.250000 # N1 N1
pair_coeff 2 2 0.065999 3.500000 # C2 C2
pair_coeff 3 3 0.065999 3.500000 # C3 C3
pair_coeff 4 4 0.065999 3.500000 # C4 C4
pair_coeff 5 5 0.065999 3.500000 # C5 C5
pair_coeff 6 6 0.065999 3.500000 # C6 C6
pair_coeff 7 7 0.065999 3.500000 # C7 C7
pair_coeff 8 8 0.065999 3.500000 # C8 C8
pair_coeff 9 9 0.140000 2.900000 # O9 O9
pair_coeff 10 10 0.065999 3.500000 # C10 C10
pair_coeff 11 11 0.065999 3.500000 # C11 C11
pair_coeff 12 12 0.030000 2.500000 # H12 H12
pair_coeff 13 13 0.030000 2.500000 # H13 H13
pair_coeff 14 14 0.030000 2.500000 # H14 H14
pair_coeff 15 15 0.030000 2.500000 # H15 H15
pair_coeff 16 16 0.030000 2.500000 # H16 H16
pair_coeff 17 17 0.030000 2.500000 # H17 H17
pair_coeff 18 18 0.030000 2.500000 # H18 H18
pair_coeff 19 19 0.030000 2.500000 # H19 H19
pair_coeff 20 20 0.030000 2.500000 # H20 H20
pair_coeff 21 21 0.030000 2.500000 # H21 H21
pair_coeff 22 22 0.030000 2.500000 # H22 H22
pair_coeff 23 23 0.030000 2.500000 # H23 H23
pair_coeff 24 24 0.030000 2.500000 # H24 H24
pair_coeff 25 25 0.030000 2.500000 # H25 H25
pair_coeff 26 26 0.030000 2.500000 # H26 H26
pair_coeff 27 27 0.030000 2.500000 # H27 H27
pair_coeff 28 28 0.030000 2.500000 # H28 H28
pair_coeff 29 29 0.030000 2.500000 # H29 H29
pair_coeff 30 30 0.030000 2.500000 # H30 H30
pair_coeff 31 31 0.030000 2.500000 # H31 H31
pair_coeff 32 32 0.030000 2.500000 # H32 H32
pair_coeff 33 33 0.030000 2.500000 # H33 H33
pair_coeff 34 34 0.250000 3.550000 # S1 S1
pair_coeff 35 35 0.170000 3.250000 # N2 N2
pair_coeff 36 36 0.170000 2.960000 # O3 O3
pair_coeff 37 37 0.170000 2.960000 # O4 O4
pair_coeff 38 38 0.060000 2.900000 # F6 F6
pair_coeff 39 39 0.060000 2.900000 # F7 F7
pair_coeff 40 40 0.060000 2.900000 # F8 F8
pair_coeff 41 41 0.250000 3.550000 # S9 S9
pair_coeff 42 42 0.170000 2.960000 # O10 O10
pair_coeff 43 43 0.060000 2.900000 # F12 F12
pair_coeff 44 44 0.060000 2.900000 # F13 F13
pair_coeff 45 45 0.060000 2.900000 # F14 F14
pair_coeff 46 46 0.170000 2.960000 # O15 O15
minimize 1.0e-4 1.0e-6 100 1000
reset_timestep 0
#fix SHAKE all shake 0.0001 20 0
neighbor 2.0 bin
neigh_modify delay 0 every 5 check yes page 1000000 one 20000
timestep 1.0
# QEQ: 2ps * 5000 = 10000 ps=10ns
fix 1 all qeq/point 10 10 1.0e-6 5000 param.qeq1
####################################################
# #
# 1. NH: NPT #
# #
####################################################
variable PBAR equal 1.0
variable TK equal 400.0
velocity all create ${TK} 12345
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
thermo_style custom step cpu etotal ke pe evdwl ecoul elong temp press vol density
thermo 5000
# monitor charge
dump TRAJ all custom 10000 eq.lammpstrj id mol type element q xu yu zu
dump_modify TRAJ element N C C C C C C C O C C H H H H H H H H H H H H H H H H H H H H H H S N O O F F F S O F F F O
# NPT 50 ps at 400 K
run 50000
unfix TPSTAT
variable TK equal 298.15
velocity all scale ${TK}
fix TPSTAT all npt temp ${TK} ${TK} 100 iso ${PBAR} ${PBAR} 1000
# NPT 100 ps at 298.15 K
run 100000
unfix TPSTAT
fix TSTAT all nvt temp ${TK} ${TK} 100
# NVT 300 ps
run 300000
write_data eq.data
write_restart eq.restart
####################################################
# #
# 3. NVT: Viscosity & MSD & RDF #
# #
####################################################
variable T equal 298.15
variable V equal vol
variable dt equal 1.0
# Maxdt=1 ns=2000000 fs, nunmber of trajectories = 10 ns/ 1 ns = 10
variable s equal 10 # sample interval
variable p equal 500000 # number of correlation time windows
variable d equal 5000000 # dump interval 2 ns
variable kB equal 1.3806504e-23 # [J/K] Boltzmann
variable atm2Pa equal 101325.0
variable A2m equal 1.0e-10
variable fs2s equal 1.0e-15
variable convert equal ${atm2Pa}*${atm2Pa}*${fs2s}*${A2m}*${A2m}*${A2m}
reset_timestep 0
variable pxy equal pxy
variable pxz equal pxz
variable pyz equal pyz
# Accumulate samples & Calculate autocorrelations for every dt
fix SS all ave/correlate $s $p $d &
v_pxy v_pxz v_pyz type auto file S0St.dat ave running
# Perform the GK running integral
variable scale equal ${convert}/(${kB}*$T)*$V*$s*${dt}
variable v11 equal trap(f_SS[3])*${scale}
variable v22 equal trap(f_SS[4])*${scale}
variable v33 equal trap(f_SS[5])*${scale}
variable v equal (v_v11+v_v22+v_v33)/3.0
variable ndens equal count(all)/vol
timestep ${dt}
fix extra1 all ave/time 10 1 10 v_pxy v_pxz v_pyz file pxyz.txt
fix extra2 all ave/time $d 1 $d v_v11 v_v22 v_v33 file v123.txt
thermo_style custom step temp press v_pxy v_pxz v_pyz v_v11 v_v22 v_v33 v_v
thermo $d
restart 2000000 meta.inbetween.restart meta.inbetween.restart
run 10000000 # 10 ns
write_data compute.data
write_restart compute.restart
print "average viscosity: $v [Pa.s] @ $T K, ${ndens} /A^3"
```

Thank you very much!