Dear Gromacs users,
I am doing protein in lipid-bilayer simulation and i am following the
procedure as per justin tutorial. I am able to insert the protein in
lipid bilayer and minimize the system as per Inflategro
procedure,during the total procedure the system was minimized in every
step.Then, I solvated and ionized sytem and minimized using the
following mdp file:
; minim.mdp - used as input into grompp to generate em.tpr
; Parameters describing what to do, when to stop and what to save
define = -DSTRONG_POSRES
integrator = steep ; Algorithm (steep = steepest descent
minimization)
emtol = 500.0 ; Stop minimization when the maximum
force < 1000.0 kJ/mol/nm
emstep = 0.01 ; Energy step size
nsteps = 50000 ; Maximum number of (minimization)
steps to perform
; Parameters describing how to find the neighbors of each atom and how
to calculate the interactions
nstlist = 1 ; Frequency to update the neighbor
list and long range forces
ns_type = grid ; Method to determine neighbor list
(simple, grid)
rlist = 1.2 ; Cut-off for making neighbor list
(short range forces)
coulombtype = PME ; Treatment of long range
electrostatic interactions
rcoulomb = 1.2 ; Short-range electrostatic cut-off
rvdw = 1.2 ; Short-range Van der Waals cut-off
pbc = xyz ; Periodic Boundary Conditions (yes/no)
the output was as follows:
Steepest Descents converged to Fmax < 500 in 4770 steps
Potential Energy = -3.8820288e+05
Maximum force = 4.4803549e+02 on atom 3573
Norm of force = 1.7854408e+01
As the potential energy and Fmax values were agreeable , I proceeded
to equillibrate the system using NVT.
The mdp file used for NVT equillibration is :
title = NVT equilibration
define = -DPOSRES ; position restrain the protein
; Run parameters
integrator = md ; leap-frog integrator
nsteps = 50000 ; 2 * 50000 = 100 ps
dt = 0.002 ; 2 fs
; Output control
nstxout = 100 ; save coordinates every 0.2 ps
nstvout = 100 ; save velocities every 0.2 ps
nstenergy = 100 ; save energies every 0.2 ps
nstlog = 100 ; update log file every 0.2 ps
; Bond parameters
continuation = no ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = all-bonds ; all bonds (even heavy atom-H bonds)
constrained
lincs_iter = 1 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
; Neighborsearching
ns_type = grid ; search neighboring grid cels
nstlist = 5 ; 10 fs
rlist = 1.2 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.2 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range
electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.16 ; grid spacing for FFT
; Temperature coupling is on
tcoupl = V-rescale ; modified Berendsen thermostat
tc-grps = Protein DPPC SOL_CL- ; three coupling groups - more accurate
tau_t = 0.1 0.1 0.1 ; time constant, in ps
ref_t = 323 323 323 ; reference temperature, one
for each group, in K
; Pressure coupling is off
pcoupl = no ; no pressure coupling in NVT
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes ; assign velocities from Maxwell distribution
gen_temp = 323 ; temperature for Maxwell distribution
gen_seed = -1 ; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 1
comm-mode = Linear
comm-grps = Protein_DPPC SOL_CL-
and the output error for the mdrun is as under:
Step 0, time 0 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms 0.002307, max 0.080808 (between atoms 3569 and 3570)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
starting mdrun 'Protein in DPPC in water'
50000 steps, 100.0 ps.
Step 0, time 0 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms 31.024508, max 1431.875854 (between atoms 448 and 450)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
453 455 49.9 0.1331 0.2004 0.1330
453 454 50.5 0.1231 0.1762 0.1230
451 452 89.9 0.1003 29.4775 0.1000
450 451 88.8 0.1364 173.8469 0.1360
448 450 89.4 0.1395 199.1697 0.1390
448 449 87.6 0.1093 41.7673 0.1090
446 450 89.5 0.1396 171.8707 0.1390
446 447 86.7 0.1093 29.7511 0.1090
444 448 84.1 0.1393 29.2753 0.1390
444 445 85.6 0.1093 16.5969 0.1090
442 446 83.9 0.1393 31.1317 0.1390
442 443 81.1 0.1093 7.2542 0.1090
441 444 84.4 0.1396 16.9881 0.1390
441 442 75.4 0.1396 5.1277 0.1390
440 441 54.3 0.1534 4.4912 0.1530
439 453 71.5 0.1532 0.8628 0.1530
439 440 85.1 0.1535 3.5898 0.1530
437 439 73.6 0.1472 0.8916 0.1470
437 438 55.1 0.1001 0.1549 0.1000
435 437 50.7 0.1331 0.2031 0.1330
3589 3590 34.2 0.1232 0.1707 0.1230
3588 3589 148.5 0.1365 0.9604 0.1360
3577 3578 40.2 0.1532 0.2055 0.1530
3576 3577 162.3 0.1532 0.1352 0.1530
3575 3576 60.7 0.1536 2.3589 0.1530
3574 3575 127.4 0.1602 4.5642 0.1530
3573 3574 72.3 0.1583 43.2879 0.1530
3572 3573 95.7 0.1638 186.7845 0.1530
3570 3572 82.5 0.1561 179.5817 0.1530
3570 3571 110.9 0.1297 27.9543 0.1230
3569 3570 127.4 0.1470 29.7186 0.1360
3568 3587 143.4 0.1552 8.5459 0.1530
3567 3568 107.8 0.1542 8.2206 0.1530
3566 3567 37.6 0.1432 0.4594 0.1430
3563 3566 115.0 0.1610 0.1509 0.1610
step 0Warning: 1-4 interaction between 435 and 440 at distance 3.609
which is larger than the 1-4 table size 2.200 nm
These are ignored for the rest of the simulation
This usually means your system is exploding,
if not, you should increase table-extension in your mdp file
or with user tables increase the table size
Segmentation fault
I have used a minimized protein as a starting structure for insertion,
and the coordinates of the atoms in the pdb file are as follows:
ATOM 435 CD ILE 45 59.176 46.081 35.828 1.00 0.00
ATOM 436 C ILE 45 55.567 47.340 36.239 1.00 0.00
ATOM 437 O ILE 45 55.953 47.123 37.392 1.00 0.00
ATOM 438 N TYR 46 54.494 48.075 35.978 1.00 0.00
ATOM 439 H TYR 46 54.293 48.397 35.053 1.00 0.00
ATOM 440 CA TYR 46 53.607 48.618 37.037 1.00 0.00
and my nvt.log file is as foolows;
:-) mdrun (-:
Input Parameters:
integrator = md
nsteps = 50000
init_step = 0
ns_type = Grid
nstlist = 5
ndelta = 2
nstcomm = 1
comm_mode = Linear
nstlog = 100
nstxout = 100
nstvout = 100
nstfout = 0
nstenergy = 100
nstxtcout = 0
init_t = 0
delta_t = 0.002
xtcprec = 1000
nkx = 48
nky = 48
nkz = 60
pme_order = 4
ewald_rtol = 1e-05
ewald_geometry = 0
epsilon_surface = 0
optimize_fft = FALSE
ePBC = xyz
bPeriodicMols = FALSE
bContinuation = FALSE
bShakeSOR = FALSE
etc = V-rescale
epc = No
epctype = Isotropic
tau_p = 1
ref_p (3x3):
ref_p[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
ref_p[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
ref_p[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
compress (3x3):
compress[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
compress[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
compress[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
refcoord_scaling = No
posres_com (3):
posres_com[0]= 0.00000e+00
posres_com[1]= 0.00000e+00
posres_com[2]= 0.00000e+00
posres_comB (3):
posres_comB[0]= 0.00000e+00
posres_comB[1]= 0.00000e+00
posres_comB[2]= 0.00000e+00
andersen_seed = 815131
rlist = 1.2
rtpi = 0.05
coulombtype = PME
rcoulomb_switch = 0
rcoulomb = 1.2
vdwtype = Cut-off
rvdw_switch = 0
rvdw = 1.2
epsilon_r = 1
epsilon_rf = 1
tabext = 1
implicit_solvent = No
gb_algorithm = Still
gb_epsilon_solvent = 80
nstgbradii = 1
rgbradii = 2
gb_saltconc = 0
gb_obc_alpha = 1
gb_obc_beta = 0.8
gb_obc_gamma = 4.85
sa_surface_tension = 2.092
DispCorr = EnerPres
free_energy = no
init_lambda = 0
sc_alpha = 0
sc_power = 0
sc_sigma = 0.3
delta_lambda = 0
nwall = 0
wall_type = 9-3
wall_atomtype[0] = -1
wall_atomtype[1] = -1
wall_density[0] = 0
wall_density[1] = 0
wall_ewald_zfac = 3
pull = no
disre = No
disre_weighting = Conservative
disre_mixed = FALSE
dr_fc = 1000
dr_tau = 0
nstdisreout = 100
orires_fc = 0
orires_tau = 0
nstorireout = 100
dihre-fc = 1000
em_stepsize = 0.01
em_tol = 10
niter = 20
fc_stepsize = 0
nstcgsteep = 1000
nbfgscorr = 10
ConstAlg = Lincs
shake_tol = 1e-04
lincs_order = 4
lincs_warnangle = 30
lincs_iter = 1
bd_fric = 0
ld_seed = 1993
cos_accel = 0
deform (3x3):
deform[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
deform[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
deform[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
userint1 = 0
userint2 = 0
userint3 = 0
userint4 = 0
userreal1 = 0
userreal2 = 0
userreal3 = 0
userreal4 = 0
grpopts:
nrdf: 6048.01 12219 45018
ref_t: 323 323 323
tau_t: 0.1 0.1 0.1
anneal: No No No
ann_npoints: 0 0 0
acc: 0 0 0
nfreeze: N N N
energygrp_flags[ 0]: 0
efield-x:
n = 0
efield-xt:
n = 0
efield-y:
n = 0
efield-yt:
n = 0
efield-z:
n = 0
efield-zt:
n = 0
bQMMM = FALSE
QMconstraints = 0
QMMMscheme = 0
scalefactor = 1
qm_opts:
ngQM = 0
Table routines are used for coulomb: TRUE
Table routines are used for vdw: FALSE
Will do PME sum in reciprocal space.
++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
U. Essman, L. Perela, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen
A smooth particle mesh Ewald method
J. Chem. Phys. 103 (1995) pp. 8577-8592
-------- -------- --- Thank You --- -------- --------
Using a Gaussian width (1/beta) of 0.384195 nm for Ewald
Cut-off's: NS: 1.2 Coulomb: 1.2 LJ: 1.2
System total charge: 0.000
Generated table with 1100 data points for Ewald.
Tabscale = 500 points/nm
Generated table with 1100 data points for LJ6.
Tabscale = 500 points/nm
Generated table with 1100 data points for LJ12.
Tabscale = 500 points/nm
Generated table with 1100 data points for 1-4 COUL.
Tabscale = 500 points/nm
Generated table with 1100 data points for 1-4 LJ6.
Tabscale = 500 points/nm
Generated table with 1100 data points for 1-4 LJ12.
Tabscale = 500 points/nm
Enabling SPC water optimization for 7497 molecules.
Configuring nonbonded kernels...
Testing AMD 3DNow support... not present.
Testing ia32 SSE support... present.
Removing pbc first time
Initializing LINear Constraint Solver
++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
B. Hess and H. Bekker and H. J. C. Berendsen and J. G. E. M. Fraaije
LINCS: A Linear Constraint Solver for molecular simulations
J. Comp. Chem. 18 (1997) pp. 1463-1472
-------- -------- --- Thank You --- -------- --------
The number of constraints is 9045
++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
S. Miyamoto and P. A. Kollman
SETTLE: An Analytical Version of the SHAKE and RATTLE Algorithms for Rigid
Water Models
J. Comp. Chem. 13 (1992) pp. 952-962
-------- -------- --- Thank You --- -------- --------
Center of mass motion removal mode is Linear
We have the following groups for center of mass motion removal:
0: Protein_DPPC
1: SOL_CL-
++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
G. Bussi, D. Donadio and M. Parrinello
Canonical sampling through velocity rescaling
J. Chem. Phys. 126 (2007) pp. 014101
-------- -------- --- Thank You --- -------- --------
There are: 31609 Atoms
Max number of connections per atom is 28
Total number of connections is 138849
Max number of graph edges per atom is 4
Total number of graph edges is 48078
Constraining the starting coordinates (step 0)
Constraining the coordinates at t0-dt (step 0)
RMS relative constraint deviation after constraining: 6.80e-04
Initial temperature: 338.401 K
Started mdrun on node 0 Sun Jul 5 17:09:45 2009
Step Time Lambda
0 0.00000 0.00000
Grid: 10 x 10 x 13 cells
Long Range LJ corr.: <C6> 9.7327e-04
Long Range LJ corr.: Epot -2200.19, Pres: -136.406, Vir: 2200.19
Energies (kJ/mol)
Angle G96Angle Proper Dih. Ryckaert-Bell. Improper Dih.
2.20077e+04 8.54042e+03 6.78950e+03 4.34650e+03 3.03266e+03
LJ-14 Coulomb-14 LJ (SR) Disper. corr. Coulomb (SR)
4.76527e+03 5.50236e+04 8.36617e+06 -2.20019e+03 -4.46844e+05
Coul. recip. Position Rest. Potential Kinetic En. Total Energy
-1.65524e+05 1.07769e+01 7.85612e+06 5.88813e+10 5.88892e+10
Conserved En. Temperature Pressure (bar) Cons. rmsd ()
5.88892e+10 2.23805e+08 1.21713e+09 3.10245e+01
Please help me to proceed further and let me know where are the
mistakes lying and how to overcome them.
Thanks in advance,
Ram
_______________________________________________
gmx-users mailing list gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
