The simulation results of a freely jointed chain do not match the theoretical values

I am new to LAMMPS and preparing to work on polymer simulations, so I chose a small system for practice.

I simulated freely jointed chains with lengths of 50, 100, 200, and 400. I only set bond_style to harmonic, while all other styles for pair, angle, and dihedral were set to none, ensuring that only bond length constraints affect the free chain.

The theoretical mean square end-to-end distance is given by R_g^2 = n * l^2. However, apart from the chain with a length of 50, the simulated end-to-end distances deviate significantly from the theoretical values. I am unsure what the issue is.

I have already tried adjusting the timestep, temperature, and ensemble, but none of these changes improved the results. I also analyzed the bond lengths and bond angles from the lammpstrj file, and they follow the theoretical distribution, ruling out any issues related to randomness.

Attached are the in and data files for the chain length of 100, as well as the data file for the chain length of 200. The theoretical end-to-end distances for a chain length of 4.48 are 44.8 and 63.36, but the actual results are only 33 and 39.
in.PDMS
PDMS_100.data
PDMS_200.data

It is not obvious to me that a 100-mer chain can fully equilibrate in only one million time steps.

Following @srtee’s very valid comment, I would suggest the following method to assess whether your system has equilibrated: measure the end-to-end distance over time and ensure that it has converged.

I don’t know your system, but for some polymer systems, it can take the microsecond or more.