[gmx-users] free energy perturbation segmentation fault.
Francesca Vitalini
francesca.vitalini11 at gmail.com
Fri Jan 24 11:29:41 CET 2014
Dear GROMACS users,
I am trying to perform a perturbation of a valine Ac-VAL-NHMe into an
alanine Ac-ALA-NHMe. Recently I have posted on the gromacs mailing
list [Free energy perturbation parameters sc-alpha sc-power sc-sigma ]
as I was experiencing problems (LINCS warning and crashes) in equilibrating
the system at lambda=.08 and lambda=1.
I am glad to say that using as a template the mdp files from the Free
energy calculations Tutorial by Justin Lemkul and reducing the timestep I
have managed to equilibrate successfully (NVT and NPT) the case lambda=0.8.
However, when trying to equilibrate at lambda=1 I end up with a
"Segmentation Fault" error which I am not able to interpret.
Does anyone have suggestions about it?
Here I summarize my problem again. I am trying to simulate the
transformation Valine to Alanine using 6 different lambda points
(0,0.2,0.4,0.6,0.8,1). At each of them I would like to have an equilibrated
system to simulate at fixed lambda.
Thank you very much.
Francesca
*FILES*:
* Here is my topology, where I make the HG of valine disappear and the CG
transform into HB.
; Include forcefield parameters
#include "oplsaa.ff/forcefield.itp"
[ moleculetype ]
; Name nrexcl
Protein_chain_A 3
[ atoms ]
; nr type resnr residue atom cgnr charge mass
typeB chargeB massB
; residue 1 ACE rtp ACE q 0.0
1 opls_135 1 ACE CH3 1 -0.18 12.011
opls_135 -0.18 12.011
2 opls_140 1 ACE HH31 1 0.06 1.008
opls_140 0.06 1.008
3 opls_140 1 ACE HH32 1 0.06 1.008
opls_140 0.06 1.008
4 opls_140 1 ACE HH33 1 0.06 1.008
opls_140 0.06 1.008
5 opls_235 1 ACE C 2 0.5 12.011
opls_235 0.5 12.011
6 opls_236 1 ACE O 2 -0.5 15.9994
opls_236 -0.5 15.9994
; residue 2 VAL rtp VAL q 0.0
7 opls_238 2 VAL N 3 -0.5 14.0067
opls_238 -0.5 14.0067
8 opls_241 2 VAL H 3 0.3 1.008
opls_241 0.3 1.008
9 opls_224B 2 VAL CA 3 0.14 12.011
opls_224B 0.14 12.011
10 opls_140 2 VAL HA 3 0.06 1.008
opls_140 0.06 1.008
11 opls_137 2 VAL CB 4 -0.06 12.011
opls_137 -0.18 12.011
12 opls_140 2 VAL HB 4 0.06 1.008
opls_140 0.06 1.008
13 opls_135 2 VAL CG1 5 -0.18 12.011
opls_140 0.06 1.008
14 opls_140 2 VAL HG11 5 0.06 1.008
opls_140 0.00 0.000
15 opls_140 2 VAL HG12 5 0.06 1.008
opls_140 0.00 0.000
16 opls_140 2 VAL HG13 5 0.06 1.008
opls_140 0.00 0.000
17 opls_135 2 VAL CG2 6 -0.18 12.011
opls_140 0.06 1.008
18 opls_140 2 VAL HG21 6 0.06 1.008
opls_140 0.00 0.000
19 opls_140 2 VAL HG22 6 0.06 1.008
opls_140 0.00 0.000
20 opls_140 2 VAL HG23 6 0.06 1.008
opls_140 0.00 0.000
21 opls_235 2 VAL C 7 0.5 12.011
opls_235 0.5 12.011
22 opls_236 2 VAL O 7 -0.5 15.9994
opls_236 -0.5 15.9994
; residue 3 NAC rtp NAC q 0.0
23 opls_238 3 NAC N 8 -0.5 14.0067
opls_238 -0.5 14.0067
24 opls_241 3 NAC H 8 0.3 1.008
opls_241 0.3 1.008
25 opls_242 3 NAC CH3 9 0.02 12.011
opls_242 0.02 12.011
26 opls_140 3 NAC HH31 9 0.06 1.008
opls_140 0.06 1.008
27 opls_140 3 NAC HH32 9 0.06 1.008
opls_140 0.06 1.008
28 opls_140 3 NAC HH33 9 0.06 1.008
opls_140 0.06 1.008
[ bonds ]
; ai aj funct c0 c1 c2 c3
1 2 1
1 3 1
1 4 1
1 5 1
5 6 1
5 7 1
7 8 1
7 9 1
9 10 1
9 11 1
9 21 1
11 12 1
11 13 1
11 17 1
13 14 1
13 15 1
13 16 1
17 18 1
17 19 1
17 20 1
21 22 1
21 23 1
23 24 1
23 25 1
25 26 1
25 27 1
25 28 1
[ pairs ]
; ai aj funct c0 c1 c2 c3
1 8 1
1 9 1
2 6 1
2 7 1
3 6 1
3 7 1
4 6 1
4 7 1
5 10 1
5 11 1
5 21 1
6 8 1
6 9 1
7 12 1
7 13 1
7 17 1
7 22 1
7 23 1
8 10 1
8 11 1
8 21 1
9 14 1
9 15 1
9 16 1
9 18 1
9 19 1
9 20 1
9 24 1
9 25 1
10 12 1
10 13 1
10 17 1
10 22 1
10 23 1
11 22 1
11 23 1
12 14 1
12 15 1
12 16 1
12 18 1
12 19 1
12 20 1
12 21 1
13 18 1
13 19 1
13 20 1
13 21 1
14 17 1
15 17 1
16 17 1
17 21 1
21 26 1
21 27 1
21 28 1
22 24 1
22 25 1
24 26 1
24 27 1
24 28 1
[ angles ]
; ai aj ak funct c0 c1
c2 c3
2 1 3 1
2 1 4 1
2 1 5 1
3 1 4 1
3 1 5 1
4 1 5 1
1 5 6 1
1 5 7 1
6 5 7 1
5 7 8 1
5 7 9 1
8 7 9 1
7 9 10 1
7 9 11 1
7 9 21 1
10 9 11 1
10 9 21 1
11 9 21 1
9 11 12 1
9 11 13 1
9 11 17 1
12 11 13 1
12 11 17 1
13 11 17 1
11 13 14 1
11 13 15 1
11 13 16 1
14 13 15 1
14 13 16 1
15 13 16 1
11 17 18 1
11 17 19 1
11 17 20 1
18 17 19 1
18 17 20 1
19 17 20 1
9 21 22 1
9 21 23 1
22 21 23 1
21 23 24 1
21 23 25 1
24 23 25 1
23 25 26 1
23 25 27 1
23 25 28 1
26 25 27 1
26 25 28 1
27 25 28 1
[ dihedrals ]
; ai aj ak al funct c0 c1
c2 c3 c4 c5
2 1 5 6 3
2 1 5 7 3
3 1 5 6 3
3 1 5 7 3
4 1 5 6 3
4 1 5 7 3
1 5 7 8 3
1 5 7 9 3
6 5 7 8 3
6 5 7 9 3
5 7 9 10 3
5 7 9 11 3
5 7 9 21 3
8 7 9 10 3
8 7 9 11 3
8 7 9 21 3
7 9 11 13 3 dih_VAL_chi1_N_C_C_C
7 9 11 17 3 dih_VAL_chi1_N_C_C_C
21 9 11 13 3 dih_VAL_chi1_C_C_C_CO
21 9 11 17 3 dih_VAL_chi1_C_C_C_CO
7 9 11 12 3
10 9 11 12 3
10 9 11 13 3
10 9 11 17 3
21 9 11 12 3
7 9 21 22 3
7 9 21 23 3
10 9 21 22 3
10 9 21 23 3
11 9 21 22 3
11 9 21 23 3
9 11 13 14 3
9 11 13 15 3
9 11 13 16 3
12 11 13 14 3
12 11 13 15 3
12 11 13 16 3
17 11 13 14 3
17 11 13 15 3
17 11 13 16 3
9 11 17 18 3
9 11 17 19 3
9 11 17 20 3
12 11 17 18 3
12 11 17 19 3
12 11 17 20 3
13 11 17 18 3
13 11 17 19 3
13 11 17 20 3
9 21 23 24 3
9 21 23 25 3
22 21 23 24 3
22 21 23 25 3
21 23 25 26 3
21 23 25 27 3
21 23 25 28 3
24 23 25 26 3
24 23 25 27 3
24 23 25 28 3
[ dihedrals ]
; ai aj ak al funct c0 c1
c2 c3
1 7 5 6 1 improper_O_C_X_Y
5 9 7 8 1 improper_Z_N_X_Y
9 23 21 22 1 improper_O_C_X_Y
21 25 23 24 1 improper_Z_N_X_Y
; Include Position restraint file
#ifdef POSRES
#include "Ac_V_NHMe/Ac_V_NHMe.itp"
#endif
; Include water topology
#include "oplsaa.ff/tip3p.itp"
#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
; i funct fcx fcy fcz
1 1 1000 1000 1000
#endif
; Include topology for ions
#include "oplsaa.ff/ions.itp"
[ system ]
; Name
ACETYL-VALINE-METHYLAMIDE in water
[ molecules ]
; Compound #mols
Protein_chain_A 1
SOL 682
*These instead are the mdp files I am using
NVT_VDW-Q.MDP
; Run control
integrator = sd ; Langevin dynamics
tinit = 0
dt = 0.0001
nsteps = 200000 ; 20 ps
nstcomm = 100
; Output control
nstxout = 200000
nstvout = 200000
nstfout = 0
nstlog = 20000
nstenergy = 10000
nstxtcout = 10000
xtc-precision = 10000
xtc_grps = Protein
; Neighborsearching and short-range nonbonded interactions
nstlist = 5
ns_type = grid
pbc = xyz
rlist = 1.0
; Electrostatics
coulombtype = PME
rcoulomb = 1.0
; van der Waals
vdw-type = switch
rvdw-switch = 0.8
rvdw = 0.9
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = EnerPres
; Spacing for the PME/PPPM FFT grid
fourierspacing = 0.12
; EWALD/PME/PPPM parameters
pme_order = 6
ewald_rtol = 1e-06
epsilon_surface = 0
optimize_fft = no
; Temperature coupling
; tcoupl is implicitly handled by the sd integrator
tc_grps = system
tau_t = 1.0
ref_t = 300
; Pressure coupling is off for NVT
Pcoupl = No
tau_p = 0.5
compressibility = 4.5e-05
ref_p = 1.0
; Free energy control stuff
free_energy = yes
init_lambda = 1
delta_lambda = 0
;foreign_lambda = 0.2
sc-alpha = 0
sc-sigma = 0
couple-moltype = Protein_chain_A ; name of moleculetype to
decouple
couple-lambda0 = vdw-q ; only van der Waals interactions
couple-lambda1 = vdw ; turn off everything, in this case only
vdW
couple-intramol = no
nstdhdl = 10
; Generate velocities to start
gen_vel = yes
gen_temp = 300
gen_seed = -1
; options for bonds
constraints = h-bonds ; we only have C-H bonds here
; Type of constraint algorithm
constraint-algorithm = lincs
; Do not constrain the starting configuration
continuation = no
; Highest order in the expansion of the constraint coupling matrix
lincs-order = 12
NVT_VDW.MDP [only Free energy parts. The rest is the same as before]
; Free energy control stuff
free_energy = yes
init_lambda = 1
delta_lambda = 0
;foreign_lambda = 0.2
sc-alpha = 0.5
sc-power = 1.0
sc-sigma = 0.3
couple-moltype = Protein_chain_A ; name of moleculetype to
decouple
couple-lambda0 = vdw ; only van der Waals interactions
couple-lambda1 = none ; turn off everything, in this case
only vdW
couple-intramol = no
*This instead is the log file coming from:
grompp -f nvt_vdw-q.mdp -c em_l-bfgs.gro -p V2A.top -o nvt_vdw-q.tpr
-maxwarn 37
Linking all bonded interactions to atoms
There are 344 inter charge-group exclusions,
will use an extra communication step for exclusion forces for PME
The initial number of communication pulses is: X 2 Y 1
The initial domain decomposition cell size is: X 0.69 nm Y 1.39 nm
The maximum allowed distance for charge groups involved in interactions is:
non-bonded interactions 1.000 nm
two-body bonded interactions (-rdd) 1.000 nm
multi-body bonded interactions (-rdd) 0.694 nm
When dynamic load balancing gets turned on, these settings will change to:
The maximum number of communication pulses is: X 2 Y 1
The minimum size for domain decomposition cells is 0.555 nm
The requested allowed shrink of DD cells (option -dds) is: 0.80
The allowed shrink of domain decomposition cells is: X 0.80 Y 0.72
The maximum allowed distance for charge groups involved in interactions is:
non-bonded interactions 1.000 nm
two-body bonded interactions (-rdd) 1.000 nm
multi-body bonded interactions (-rdd) 0.555 nm
*The output from grompp is
Generated 332520 of the 332520 non-bonded parameter combinations
Generating 1-4 interactions: fudge = 0.5
Generated 332520 of the 332520 1-4 parameter combinations
WARNING 1 [file V2A.top, line 59]:
No default Bond types for perturbed atoms, using normal values
WARNING 2 [file V2A.top, line 60]:
No default Bond types for perturbed atoms, using normal values
WARNING 3 [file V2A.top, line 61]:
No default Bond types for perturbed atoms, using normal values
WARNING 4 [file V2A.top, line 62]:
No default Bond types for perturbed atoms, using normal values
WARNING 5 [file V2A.top, line 63]:
No default Bond types for perturbed atoms, using normal values
WARNING 6 [file V2A.top, line 64]:
No default Bond types for perturbed atoms, using normal values
WARNING 7 [file V2A.top, line 161]:
No default Angle types for perturbed atoms, using normal values
WARNING 8 [file V2A.top, line 162]:
No default Angle types for perturbed atoms, using normal values
WARNING 9 [file V2A.top, line 163]:
No default Angle types for perturbed atoms, using normal values
WARNING 10 [file V2A.top, line 164]:
No default Angle types for perturbed atoms, using normal values
WARNING 11 [file V2A.top, line 165]:
No default Angle types for perturbed atoms, using normal values
WARNING 12 [file V2A.top, line 166]:
No default Angle types for perturbed atoms, using normal values
WARNING 13 [file V2A.top, line 167]:
No default Angle types for perturbed atoms, using normal values
WARNING 14 [file V2A.top, line 168]:
No default Angle types for perturbed atoms, using normal values
WARNING 15 [file V2A.top, line 169]:
No default Angle types for perturbed atoms, using normal values
WARNING 16 [file V2A.top, line 170]:
No default Angle types for perturbed atoms, using normal values
WARNING 17 [file V2A.top, line 171]:
No default Angle types for perturbed atoms, using normal values
WARNING 18 [file V2A.top, line 172]:
No default Angle types for perturbed atoms, using normal values
WARNING 19 [file V2A.top, line 204]:
Some parameters for bonded interaction involving perturbed atoms are
specified explicitly in state A, but not B - copying A to B
WARNING 20 [file V2A.top, line 219]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 21 [file V2A.top, line 220]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 22 [file V2A.top, line 221]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 23 [file V2A.top, line 222]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 24 [file V2A.top, line 223]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 25 [file V2A.top, line 224]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 26 [file V2A.top, line 225]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 27 [file V2A.top, line 226]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 28 [file V2A.top, line 227]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 29 [file V2A.top, line 228]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 30 [file V2A.top, line 229]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 31 [file V2A.top, line 230]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 32 [file V2A.top, line 231]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 33 [file V2A.top, line 232]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 34 [file V2A.top, line 233]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 35 [file V2A.top, line 234]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 36 [file V2A.top, line 235]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
WARNING 37 [file V2A.top, line 236]:
No default Ryckaert-Bell. types for perturbed atoms, using normal values
Excluding 3 bonded neighbours molecule type 'Protein_chain_A'
turning H bonds into constraints...
Excluding 2 bonded neighbours molecule type 'SOL'
turning H bonds into constraints...
Coupling 1 copies of molecule type 'Protein_chain_A'
Setting gen_seed to 29390
Velocities were taken from a Maxwell distribution at 300 K
Analysing residue names:
There are: 3 Protein residues
There are: 682 Water residues
Analysing Protein...
Number of degrees of freedom in T-Coupling group System is 4157.00
Largest charge group radii for Van der Waals: 0.148, 0.103 nm
Largest charge group radii for Coulomb: 0.148, 0.148 nm
NOTE 1 [file nvt_vdw-q.mdp]:
The sum of the two largest charge group radii (0.251266) is larger than
rlist (1.000000) - rvdw (0.900000)
Calculating fourier grid dimensions for X Y Z
Using a fourier grid of 24x24x24, spacing 0.116 0.116 0.116
Estimate for the relative computational load of the PME mesh part: 0.69
NOTE 2 [file nvt_vdw-q.mdp]:
The optimal PME mesh load for parallel simulations is below 0.5
and for highly parallel simulations between 0.25 and 0.33,
for higher performance, increase the cut-off and the PME grid spacing
This run will generate roughly 1 Mb of data
There were 2 notes
There were 37 warnings
More information about the gromacs.org_gmx-users
mailing list