Dear Gromacs users I am doing free energy calculation of the protein-ligand system in gromacs 2020 in which I am annihilating ligand by removing charges and vdw interactions. During Vdw removal I am facing lincs warning and my system crashes by generating different pdb structures. I have tried to use constraints on all bonds and also reduced time step to 1fs but still getting same problem as below:-
Overriding thread affinity set outside gmx mdrun WARNING: There are no atom pairs for dispersion correction starting mdrun 'GROtesk MACabre and Sinister in water' 1000000 steps, 1000.0 ps. Step 853159, time 853.159 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.000569, max 0.041969 (between atoms 3430 and 3439) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 3435 3437 38.0 0.1400 0.1415 0.1400 3437 3438 52.2 0.1080 0.1068 0.1080 3437 3439 39.0 0.1400 0.1402 0.1400 Step 853159, time 853.159 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.000550, max 0.040414 (between atoms 3430 and 3439) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 3435 3437 37.8 0.1400 0.1415 0.1400 3437 3438 51.9 0.1080 0.1070 0.1080 3437 3439 38.6 0.1400 0.1403 0.1400 Step 853160, time 853.16 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.000833, max 0.054713 (between atoms 3437 and 3438) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 3430 3439 30.5 0.1457 0.1428 0.1400 Please find the mdp file as below:- ; Run control integrator = sd ; stochastic leap-frog integrator dt = 0.001 nsteps = 1000000 ; 1 ns comm-mode = Linear ; remove center of mass translation nstcomm = 100 ; frequency for center of mass motion removal ; Output control nstxout = 5000 nstvout = 5000 nstfout = 0 nstlog = 5000 nstenergy = 5000 nstxout-compressed = 5000 ; Bond parameters continuation = yes ; formerly known as 'unconstrained-start' - useful for exact continuations and reruns constraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; bonds to H are constrained lincs_iter = 1 ; accuracy of LINCS lincs-order = 6 ; Highest order in the expansion of the constraint coupling matrix ; Neighborsearching and short-range nonbonded interactions cutoff-scheme = verlet nstlist = 10 ns_type = grid rlist = 1.2 pbc = xyz ; 3-D PBC ; Periodic boundary conditions ; Electrostatics coulombtype = PME rcoulomb = 1.2 pme-order = 6 ; interpolation order for PME (default is 4) fourierspacing = 0.10 ; grid spacing for FFT ewald-rtol = 1e-6 ; relative strength of the Ewald-shifted direct potential at rcoulomb ewald_geometry = 3d ; Ewald sum is performed in all three dimensions ; van der Waals vdwtype = cutoff vdw-modifier = Potential-shift-Verlet verlet-buffer-tolerance = 0.005 rvdw = 1.2 rvdw-switch = 1.0 DispCorr = EnerPres ; apply long range dispersion corrections for Energy and Pressure ; Temperature coupling ; tcoupl is implicitly handled by the sd integrator tc-grps = system ; two coupling groups - more accurate tau_t = 0.1 ; time constant, in ps ref_t = 298.15 ; Pressure coupling is on for NPT Pcoupl = Parrinello-Rahman pcoupltype = isotropic ; uniform scaling of box vectors tau_p = 2.0 ; time constant (ps) ref_p = 1.0 ; reference pressure (bar) compressibility = 4.5e-05 ; isothermal compressibility of water (bar^-1) ; velocities gen_vel = no ; Velocity generation is on (if gen_vel is 'yes', continuation should be 'no') gen_seed = -1 ; Use random seed gen_temp = 298.15 ; Free energy control stuff free_energy = yes init_lambda_state = 17 delta_lambda = 0 calc_lambda_neighbors =-1 ; only immediate neighboring windows ; Vectors of lambda specified here ; Each combination is an index that is retrieved from init_lambda_state for each simulation vdw_lambdas = 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00 coul_lambdas = 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 restraint_lambdas = 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 ; Options for the decoupling sc-alpha = 0.5 sc-coul = no ; linear interpolation of Coulomb (none in this case) sc-power = 1 sc-sigma = 0.3 nstdhdl = 100 dhdl-print-energy = yes Could you please suggest me how should I solve this problem? Thanks in advance. Sadaf -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.