Error in minimization run for free energy calculation

GROMACS version:2025.2
GROMACS modification: No

I did an minimization run for my free energy calculation. Prior to that i used couple-lambda0 = vdw # with this setting everything worked. But now i want to do the same run with: couple-lambda0 = vdw-q # but encountered an error this time:

Steepest Descents:
Tolerance (Fmax) = 1.00000e+02
Number of steps = 5000
Step= 0, Dmax= 1.0e-02 nm, Epot= -4.03462e+04 Fmax= 1.81422e+03, atom= 10
Step= 2, Dmax= 5.0e-03 nm, Epot= -4.03761e+04 Fmax= 2.25691e+03, atom= 10
Step= 3, Dmax= 6.0e-03 nm, Epot= -4.03839e+04 Fmax= 3.36219e+03, atom= 16
Step= 5, Dmax= 3.6e-03 nm, Epot= -4.04115e+04 Fmax= 9.83459e+02, atom= 25
Step= 6, Dmax= 4.3e-03 nm, Epot= -4.04253e+04 Fmax= 3.61731e+03, atom= 25
Step= 7, Dmax= 5.2e-03 nm, Epot= -4.04402e+04 Fmax= 2.65308e+03, atom= 25
Step= 9, Dmax= 3.1e-03 nm, Epot= -4.04508e+04 Fmax= 1.04395e+03, atom= 25
Step= 10, Dmax= 3.7e-03 nm, Epot= -4.04587e+04 Fmax= 3.45840e+03, atom= 25
Step= 11, Dmax= 4.5e-03 nm, Epot= -4.04702e+04 Fmax= 1.88051e+03, atom= 25
Step= 13, Dmax= 2.7e-03 nm, Epot= -4.04770e+04 Fmax= 1.36128e+03, atom= 25
Step= 14, Dmax= 3.2e-03 nm, Epot= -4.04829e+04 Fmax= 2.44404e+03, atom= 25
Step= 15, Dmax= 3.9e-03 nm, Epot= -4.04894e+04 Fmax= 2.19656e+03, atom= 25
Step= 16, Dmax= 4.6e-03 nm, Epot= -4.04922e+04 Fmax= 3.33233e+03, atom= 25
Step= 17, Dmax= 5.6e-03 nm, Epot= -4.04978e+04 Fmax= 3.33252e+03, atom= 25
Step= 19, Dmax= 3.3e-03 nm, Epot= -4.05085e+04 Fmax= 6.83306e+02, atom= 25
Step= 20, Dmax= 4.0e-03 nm, Epot= -4.05161e+04 Fmax= 3.97660e+03, atom= 25
Step= 21, Dmax= 4.8e-03 nm, Epot= -4.05285e+04 Fmax= 1.78878e+03, atom= 25
Step= 23, Dmax= 2.9e-03 nm, Epot= -4.05335e+04 Fmax= 1.65437e+03, atom= 25
Step= 24, Dmax= 3.5e-03 nm, Epot= -4.05370e+04 Fmax= 2.49089e+03, atom= 25
Step= 25, Dmax= 4.2e-03 nm, Epot= -4.05418e+04 Fmax= 2.46773e+03, atom= 25
Step= 26, Dmax= 5.0e-03 nm, Epot= -4.05431e+04 Fmax= 3.50087e+03, atom= 25
Step= 27, Dmax= 6.0e-03 nm, Epot= -4.05470e+04 Fmax= 3.64037e+03, atom= 25
Step= 29, Dmax= 3.6e-03 nm, Epot= -4.05587e+04 Fmax= 6.78146e+02, atom= 25
Step= 30, Dmax= 4.3e-03 nm, Epot= -4.05632e+04 Fmax= 4.28882e+03, atom= 25
Step= 31, Dmax= 5.2e-03 nm, Epot= -4.05767e+04 Fmax= 1.90523e+03, atom= 25
Step= 33, Dmax= 3.1e-03 nm, Epot= -4.05813e+04 Fmax= 1.79261e+03, atom= 25
Step= 34, Dmax= 3.7e-03 nm, Epot= -4.05841e+04 Fmax= 2.65710e+03, atom= 25
Step= 35, Dmax= 4.5e-03 nm, Epot= -4.05885e+04 Fmax= 2.67344e+03, atom= 25
Step= 36, Dmax= 5.4e-03 nm, Epot= -4.05891e+04 Fmax= 3.72605e+03, atom= 25
Step= 37, Dmax= 6.4e-03 nm, Epot= -4.05926e+04 Fmax= 3.95423e+03, atom= 25
Step= 39, Dmax= 3.9e-03 nm, Epot= -4.06059e+04 Fmax= 6.81430e+02, atom= 25
Step= 40, Dmax= 4.6e-03 nm, Epot= -4.06099e+04 Fmax= 4.65331e+03, atom= 25
Step= 41, Dmax= 5.6e-03 nm, Epot= -4.06258e+04 Fmax= 1.97437e+03, atom= 25
Step= 43, Dmax= 3.3e-03 nm, Epot= -4.06308e+04 Fmax= 2.03574e+03, atom= 25
Step= 44, Dmax= 4.0e-03 nm, Epot= -4.06345e+04 Fmax= 2.68072e+03, atom= 25
Step= 45, Dmax= 4.8e-03 nm, Epot= -4.06389e+04 Fmax= 3.10175e+03, atom= 25
Step= 46, Dmax= 5.8e-03 nm, Epot= -4.06420e+04 Fmax= 3.68385e+03, atom= 25
Step= 47, Dmax= 6.9e-03 nm, Epot= -4.06436e+04 Fmax= 4.65000e+03, atom= 25
Step= 48, Dmax= 8.3e-03 nm, Epot= -4.06465e+04 Fmax= 5.10755e+03, atom= 25
Step= 50, Dmax= 5.0e-03 nm, Epot= -4.06689e+04 Fmax= 9.68186e+02, atom= 25
Step= 51, Dmax= 6.0e-03 nm, Epot= -4.06785e+04 Fmax= 5.10320e+03, atom= 25
Step= 52, Dmax= 7.2e-03 nm, Epot= -4.06985e+04 Fmax= 3.72655e+03, atom= 25
Step= 53, Dmax= 8.6e-03 nm, Epot= -4.06986e+04 Fmax= 5.85831e+03, atom= 25
Step= 54, Dmax= 1.0e-02 nm, Epot= -4.07125e+04 Fmax= 7.00490e+03, atom= 25
Step= 55, Dmax= 1.2e-02 nm, Epot= -4.07323e+04 Fmax= 6.77523e+03, atom= 25
Step= 57, Dmax= 7.4e-03 nm, Epot= -4.07916e+04 Fmax= 2.96944e+03, atom= 25
Step= 58, Dmax= 8.9e-03 nm, Epot= -4.09015e+04 Fmax= 9.59393e+03, atom= 26
Step= 59, Dmax= 1.1e-02 nm, Epot= -4.11581e+04 Fmax= 2.15589e+04, atom= 26
Step= 60, Dmax= 1.3e-02 nm, Epot= -4.23807e+04 Fmax= 1.49929e+05, atom= 26
Step= 61, Dmax= 1.5e-02 nm, Epot= -4.56595e+04 Fmax= 1.01834e+06, atom= 26
Step= 64, Dmax= 4.6e-03 nm, Epot= -8.56054e+04 Fmax= 5.45279e+07, atom= 26
Step= 68, Dmax= 6.9e-04 nm, Epot= -1.29059e+05 Fmax= 2.06230e+08, atom= 26
Step= 70, Dmax= 4.2e-04 nm, Epot= -1.51067e+05 Fmax= 3.20113e+08, atom= 26
Step= 72, Dmax= 2.5e-04 nm, Epot= -3.39230e+05 Fmax= 2.31300e+09, atom= 26
Step= 75, Dmax= 7.5e-05 nm, Epot= -2.21598e+06 Fmax= 1.22024e+11, atom= 26
Step= 79, Dmax= 1.1e-05 nm, Epot= -8.27067e+06 Fmax= inf, atom= 26
Segmentation fault (core dumped)

Why is the energy going through the roof now ?

Like i said, i changed nothing except the coulomb-lambda settings:

couple-moltype = PROA ; name of moleculetype to decouple
couple-lambda0 = vdw-q ; only van der Waals interactions
couple-lambda1 = none ; turn off everything, in this case only vdW
couple-intramol = no
; Vectors of lambda specified here
; Each combination is an index that is retrieved from init_lambda_state for each simulation
; init_lambda_state 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
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 = 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
; We are not transforming any bonded or restrained interactions
bonded_lambdas = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
restraint_lambdas = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
; Masses are not changing (particle identities are the same at lambda = 0 and lambda = 1)
mass_lambdas = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00

Can someone help me ? Or does this look common to someone, If so, that would help a lot !

I am answering myself again, because found the error by myself: I haven’t changed the mdp settings i took from lemkul labs.

Because of only decoupling vdw There was:

sc-coul = no

Now when also calculating the coul_lambda states,it runs different. To fix his kind of error, you need to set:

sc-coul = yes

I dont now if this is the best or only way for this kind of error. Please correct me if iam wrong

What you’re seeing here is most likely the problem of two charged particles getting too close to each other due to the LJ interactions being decoupled, which causes the system to become unstable. The sc-coul option introduces a soft-core part to the Coulomb potential, which avoids this. As I said in the other post though, you could also consider turning the Coulomb interactions off before the LJ, which would also solve this problem.

1 Like

Interesting, thanks for the explanation and the answer here aswell ! Do you think this might be the ions i added (0,15 mol) ? Otherwise i haven’t charged my solute (TRP) manually. Yes, i’m really curious how this might work out

No, this should have nothing to do with ions. Regular atoms also have a partial charge, and if the repulsive part of the LJ potential isn’t strong enough to keep them apart anymore (because it’s scaled by lambda), charges can get too close to each other, causing your forces to explode, as you can see in your log.

I understand, makes totally sense. I was thinking of real charges and wasn’t sure how the ions will be treated here.