GROMACS version: 2022.2
GROMACS modification: No
Dear Community,
I’m using GROMACS to simulate polymers and study the stress relaxation using Green-Kubo relationship. The polymer is merely bead-spring model. I compared the results between two systems, namely the bond between beads using either a harmonic spring of force constant 200kT/b^2 or LINCS constraint. If I only include bond potential ( or LINCS) without any other potential, I get the results as shown in Figure 1.
The stress produced from LINCS algorithm didn’t really match the theory but after a vertical shift it kinda overlaps with the theory curve and that produced by harmonic spring. However when I tried to simulate the system with bending and intermolecular potentials, I got completely different results on longer time scale. (Figure 2)
This quite puzzles me. My understanding is that either harmonic spinrg or LINCS constraint should only have impact on the short time scale when each single bond relaxes. Even though they don’t overlap at the short time scale they should coincide at longer time. However my results showed the opposite. So I’m wondering how LINCS calculates stress and why this affects its relaxtion a lot.
(P.S. G(t)=\frac{V}{3k_BT} \Sigma_{xy,xz,yz}\langle \sigma_i(t) \sigma_i(0) \rangle, but I normalize it by \nu k_B T = \frac{N_{\text{chain}} k_B T}{V}, so the curve in the plot are actually only virial autocorrelation \frac{G(t)}{\nu k_B T} =\frac{1}{3N_{\text{chain}} }\Sigma_{xy,xz,yz}\langle \Xi_i(t) \Xi_i(0) \rangle where \Xi is the virial. The simulation was using Brownian solver and NVT ensemble so the kinetic contribution to the stress was ignored.)
Thank you in advance!

