Forces not calcuating properly

LAMMPS
,
1

I have made a triangular lattice. And all the points are 1 unit of distance apart and I have set the bond is harmonic with rest length of 1 so without adding any other potential or force the force on each particle should be 0 and I have calculated it using a md code and it comes out to be zero. But then I run lammps for the same configuration the forces in dump file is not zero. Here is the basic input script I am using.

dimension           2
units               lj
atom_style          bond

# Read initial configuration
read_data           init_config.txt

# Define bonds (springs)
bond_style          harmonic
bond_coeff          1 200.0 1

# Define the neighbor list
neighbor            0.3 bin
neigh_modify        delay 5
# Output atom properties
dump 1 all custom 100 trajectory/atoms.lammpstrj id type x y z
dump 2 all custom 100 forces/forces id x y z fx fy fz
# Specify simulation settings
# thermo 100

# Run the simulation for the required duration
run 1000

# Write a restart file
write_data restart.txt

And here is the output of dump file.

ITEM: TIMESTEP
0
ITEM: NUMBER OF ATOMS
100
ITEM: BOX BOUNDS pp pp pp
-5.0000000000000000e+00 1.4500000000000000e+01
-5.0000000000000000e+00 1.2794200000000000e+01
-5.0000000000000000e-01 5.0000000000000000e-01
ITEM: ATOMS id x y z fx fy fz
1 0 0 0 0 0 0
2 1 0 0 0 0 0
3 2 0 0 0 0 0
4 3 0 0 0 0 0
5 4 0 0 0 0 0
6 5 0 0 0 0 0
7 6 0 0 0 0 0
8 7 0 0 0 0 0
9 8 0 0 0 0 0
10 9 0 0 0 0 0
11 0.5 0.866025 0 0 0 0
12 1.5 0.866025 0 0 0 0
13 2.5 0.866025 0 0 0 0
14 3.5 0.866025 0 0 0 0
15 4.5 0.866025 0 0 0 0
16 5.5 0.866025 0 0 0 0
17 6.5 0.866025 0 0 0 0
18 7.5 0.866025 0 0 0 0
19 8.5 0.866025 0 0 0 0
20 9.5 0.866025 0 0 0 0
21 0 1.73205 0 0 0 0
22 1 1.73205 0 0 0.00300001 0
23 2 1.73205 0 0 0.00300001 0
24 3 1.73205 0 0 0.00300001 0
25 4 1.73205 0 0 0.00300001 0
26 5 1.73205 0 0 0.00300001 0
27 6 1.73205 0 0 0.00300001 0
28 7 1.73205 0 0 0.00300001 0
29 8 1.73205 0 0 0.00300001 0
30 9 1.73205 0 0 0 0

I don’t know what I am doing wrong here.

2

The most likely explanation would be that there is an issue with your input geometry.

3

Yeah at first I thought the same thing. So what I did was take all the positions and wrote a code to calculate forces on the particle based on the same parameters given in the input file. The forces on each one came out to be zero as it should have been.

4

Sorry, but that doesn’t prove anything. If you have a mistake in the data file, this is not being tested this way.