hi guys,

for quite a long time i am trying to 'create' a forcefield for aliphatic alcohols (right now, for ethanol) to use with a water potential we employ in our simulations.

there is a problem - the water potential uses combination rule 2, so i thought i would just try to use FFGMX and OPLS and change it a bit. it didn't seem to work. i also tried the FFAMBER, that uses exactly this combination rule. unfortunately, it fails as well.

i'm getting all sorts of errors - segfaults, lincs failures, ... in most cases, the whole simulation 'freezes' (kinetic energy of order of 1e-5, for example) and nothing happens. i am also quite frequently getting close contacts between atoms (within one molecule, or between two different molecules).

it seems i am missing something in the forcefield, but i am not able to identify what is it... i would be very grateful if anyone of you could help me and shed more light in this.

if anyone has working alcohol potential (united or all atom type) with comb. rule 2, i would be more than happy, if you could share it with me...

also a 'side question' - when using FFAMBER, there are dihedrals defined with all atom names specified (e.g., CT CT CT CT), and also dihedrals with wildcards (X CT CT X). for a CT-CT-CT-CT dihedrals, both of these would be applicable - are they both used in the simulation, or does the latter (X CT CT X) override the former definition (CT CT CT CT)?

once again, i would be very glad if anyone of you could help me. i went through the manual many times, but i am becoming more and more confused...

i attach the files for a simulation that is failing for me. it's 2+2 pentanol molecules in a prismatic simulation box (i have it setup like this because we are normally doing slab simulations). all necessary files are included.

thanks. with best regards,
lubos


;
; forcefield by Lubos Vrbka for (contaminated) ice simulations
; Lubos Vrbka, 2006
; lubos (@) vrbka (.) net
;
; water parameters - nada/van der eerden
;       implementation by marcelo carignano (cari (@) purdue (.) edu)
; nacl parameters - j. chem. phys. 100, 3757
;       implementation by marcelo carignano (cari (@) purdue (.) edu)
;       used with spc/e - consistency?
; alcohol parameters - AMBER parm99 parameters
;       implementation by lubos vrbka (lubos (@) vrbka (.) net)
;
; changed gen-pairs to yes (consistently with amber FF) - should not affect
; NE6 water potential

#define _FF_ICE

[ defaults ]
; nbfunc        comb-rule       gen-pairs       FudgeLJ FudgeQQ
  1             2               yes             0.5     0.8333

; *************************************************************************

[ atomtypes ]
;name  type     mass            charge  ptype   c6              c12
; original AMBER atom types: CT, H1, HC, OH, HO
; CT sp3 carbon (charge corresponding to CH3!)
; H1 hydrogen bound to C with el. withdrawing substituent
; HC hydrogen bound to C
; HO hydroxyl hydrogen
; OH hydroxyl oxygen 
; name  btype   mass    charge  ptype   sigma           epsilon
CT      CT       12.010  0.1116 A       3.39967e-01     4.57730e-01
HC      HC       1.0080  0.0372 A       2.64953e-01     6.56888e-02
H1      H1       1.0080  0.0372 A       2.47135e-01     6.56888e-02
;HO     HO       1.0080  0.4215 A       2.47135e-01     6.56888e-02
HO      HO       1.0080  0.4215 A       0.00000e+00     0.00000e+00
OH      OH      16.0000 -0.6497 A       3.06647e-01     8.80314e-01
; NE6 water model
EP      EP       0.00000        -0.866  D       0.00000E+00     0.00000E+00
OW      OW      15.99940         0.000  A       3.11500E-01     0.714850050
HW      HW       1.00800         0.477  A       6.73000E-02     0.115419008
LP      LP       0.00000        -0.477  D       0.00000E+00     0.00000E+00

; sections [ ????types ] contain generic alcohol contaminant params
[ bondtypes ]
; i     j       func    b0      kb
HO      OH      1       0.09600 462750.4
CT      OH      1       0.14100 267776.0
CT      CT      1       0.15260 259408.0
CT      H1      1       0.10900 284512.0
CT      HC      1       0.10900 284512.0
  
[ angletypes ]
; i     j       k       func    th0     cth
CT      OH      HO      1       108.500 460.240
CT      CT      OH      1       109.500 418.400
CT      CT      CT      1       109.500 334.720
CT      CT      HC      1       109.500 418.400
CT      CT      H1      1       109.500 418.400
HC      CT      HC      1       109.500 292.880
H1      CT      H1      1       109.500 292.880
H1      CT      OH      1       109.500 418.400

[ dihedraltypes ]
; i     j       k       l       func    coefficients
;CT     CT      OH      HO      3       1.71544  0.96232  0.00000 -2.67776  
0.00000  0.00000
;CT     CT      CT      OH
;CT     CT      CT      CT      3       3.68192  3.09616 -2.09200 -3.01248  
0.00000  0.00000
;CT     CT      CT      HC      3       0.66944  2.00832  0.00000 -2.67776  
0.00000  0.00000
;CT     CT      CT      H1      3       0.66944  2.00832  0.00000 -2.67776  
0.00000  0.00000
;H1     CT      OH      HO
;HC     CT      CT      H1
;HC     CT      CT      HC      3       0.62760  1.88280  0.00000 -2.51040  
0.00000  0.0000
;HC     CT      CT      OH      3       1.04600 -1.04600  0.00000  0.00000  
0.00000  0.00000
X       CT      CT      X       3       0.65084  1.95253  0.00000 -2.60338  
0.00000  0.00000
X       CT      OH      X       3       0.69733  2.09200  0.00000 -2.78933  
0.00000  0.00000

; *************************************************************************

[ moleculetype ]
; molname       nrexcl
5OL             3

[ atoms ]
; id    at type res nr  resname at name cg nr   charge
1       HO      1       5OL     HO      1        0.420
2       OH      1       5OL     OH      1       -0.648
3       CT      1       5OL     C01     1        0.154
4       H1      1       5OL     H11     1        0.037
5       H1      1       5OL     H12     1        0.037
6       CT      1       5OL     C02     1       -0.074
7       HC      1       5OL     H21     1        0.037
8       HC      1       5OL     H22     1        0.037
9       CT      1       5OL     C03     1       -0.074
10      HC      1       5OL     H31     1        0.037
11      HC      1       5OL     H32     1        0.037
12      CT      1       5OL     C04     1       -0.074
13      HC      1       5OL     H41     1        0.037
14      HC      1       5OL     H42     1        0.037
15      CT      1       5OL     C05     1       -0.111
16      HC      1       5OL     H51     1        0.037
17      HC      1       5OL     H52     1        0.037
18      HC      1       5OL     H53     1        0.037

[ bonds ]
; values defined after the defaults section
; i     j       funct   b0      kb
1       2       1
2       3       1
3       4       1
3       5       1
3       6       1
6       7       1
6       8       1
6       9       1
9       10      1
9       11      1
9       12      1
12      13      1
12      14      1
12      15      1
15      16      1
15      17      1
15      18      1

[ angles ]
; values defined after the defaults section
; i     j       k       funct th0       cth
1       2       3       1
2       3       4       1
2       3       5       1
2       3       6       1
3       6       7       1
3       6       8       1
3       6       9       1
6       9       10      1
6       9       11      1
6       9       12      1
9       12      13      1
9       12      14      1
9       12      15      1
12      15      16      1
12      15      17      1
12      15      18      1
4       3       5       1
4       3       6       1
7       6       8       1
7       6       9       1
10      9       11      1
10      9       12      1
13      12      14      1
13      12      15      1
16      15      17      1
16      15      18      1

[ dihedrals ]
; values defined after the defaults section
; i     j       k       l       func    coefficients
1       2       3       4       3
1       2       3       5       3
1       2       3       6       3
2       3       6       7       3
2       3       6       8       3
2       3       6       9       3
3       6       9       10      3
3       6       9       11      3
3       6       9       12      3
6       9       12      13      3
6       9       12      14      3
6       9       12      15      3
9       12      15      16      3
9       12      15      17      3
9       12      15      18      3
4       3       6       7       3
4       3       6       8       3
4       3       6       9       3
5       3       6       7       3
5       3       6       8       3
5       3       6       9       3
7       6       9       10      3
7       6       9       11      3
7       6       9       12      3
8       6       9       10      3
8       6       9       11      3
8       6       9       12      3
10      9       12      13      3
10      9       12      14      3
10      9       12      15      3
11      9       12      13      3
11      9       12      14      3
11      9       12      15      3
13      12      15      16      3
13      12      15      17      3
13      12      15      18      3
14      12      15      16      3
14      12      15      17      3
14      12      15      18      3

;[ pairs ]
; i     j       funct   c6      c12
;1      4       1
;1      5       1
;1      6       1
;1      7       1
;2      5       1
;2      6       1
;2      7       1
;3      6       1
;3      7       1
;4      7       1

; *************************************************************************

[ moleculetype ]
; molname       nrexcl
NE6             1

[ atoms ]
; We use a strange order of atoms to make things go faster in GROMACS (?)
; id    at type res nr  resname at name cg nr   charge
1       HW      1       NE6     HW1     1        0.477
2       HW      1       NE6     HW2     1        0.477
3       LP      1       NE6     LP1     1       -0.044
4       LP      1       NE6     LP2     1       -0.044
5       EP      1       NE6     EP      1       -0.866
6       OW      1       NE6     OW      1        0.000

[ constraints ]
; i     j       funct   distance
1       6       1       0.0980  ; HW1 - OW
2       6       1       0.0980  ; HW1 - OW
1       2       1       0.15857 ; HW1 - HW2

[ dummies3 ]
; For the EP:
; Dummy from                    funct   a               b
5       6       1       2       1       0.199642537     0.199642537
; For the LPs:
; Dummy from                    funct   a               b               c
3       6       1       2       4       -0.437172       -0.437172        
8.022961206
4       6       1       2       4       -0.437172       -0.437172       
-8.022961206

[ exclusions ]
1       2       3       4       5       6
2       1       3       4       5       6
3       1       2       4       5       6
4       1       2       3       5       6
5       1       2       3       4       6
6       1       2       3       4       5
;
; sample top file for NE6 water with contaminant simulation
;

; include parameters and topology
#include "_ffice.itp"

[ system ]
; Name
NE6 water with contaminant

[ molecules ]
; Compound      #mols
;NE6            192
;NA             0
;CL             0
5OL             4
;2OL            10
NE6 water with contaminant
   72
    15OL     HO    1   0.300   0.299   6.499
    15OL     OH    2   0.300   0.242   6.576
    15OL    C01    3   0.300   0.322   6.692
    15OL    H11    4   0.211   0.385   6.692
    15OL    H12    5   0.389   0.385   6.692
    15OL    C02    6   0.300   0.234   6.816
    15OL    H21    7   0.389   0.171   6.816
    15OL    H22    8   0.211   0.171   6.816
    15OL    C03    9   0.300   0.321   6.941
    15OL    H31   10   0.211   0.384   6.941
    15OL    H32   11   0.389   0.384   6.942
    15OL    C04   12   0.300   0.232   7.065
    15OL    H41   13   0.389   0.169   7.065
    15OL    H42   14   0.211   0.170   7.065
    15OL    C05   15   0.300   0.319   7.190
    15OL    H51   16   0.389   0.382   7.191
    15OL    H52   17   0.300   0.255   7.278
    15OL    H53   18   0.211   0.382   7.191
    25OL     HO   19   0.900   1.199   6.499
    25OL     OH   20   0.900   1.142   6.576
    25OL    C01   21   0.900   1.222   6.692
    25OL    H11   22   0.811   1.285   6.692
    25OL    H12   23   0.989   1.285   6.692
    25OL    C02   24   0.900   1.134   6.816
    25OL    H21   25   0.989   1.071   6.816
    25OL    H22   26   0.811   1.071   6.816
    25OL    C03   27   0.900   1.221   6.941
    25OL    H31   28   0.811   1.284   6.941
    25OL    H32   29   0.989   1.284   6.942
    25OL    C04   30   0.900   1.132   7.065
    25OL    H41   31   0.989   1.069   7.065
    25OL    H42   32   0.811   1.070   7.065
    25OL    C05   33   0.900   1.219   7.190
    25OL    H51   34   0.989   1.282   7.191
    25OL    H52   35   0.900   1.155   7.278
    25OL    H53   36   0.811   1.282   7.191
    35OL     HO   37   0.300   0.299   2.801
    35OL     OH   38   0.300   0.242   2.724
    35OL    C01   39   0.300   0.322   2.608
    35OL    H11   40   0.389   0.385   2.608
    35OL    H12   41   0.211   0.385   2.608
    35OL    C02   42   0.300   0.234   2.484
    35OL    H21   43   0.211   0.171   2.484
    35OL    H22   44   0.389   0.171   2.484
    35OL    C03   45   0.300   0.321   2.359
    35OL    H31   46   0.389   0.384   2.359
    35OL    H32   47   0.211   0.384   2.358
    35OL    C04   48   0.300   0.232   2.235
    35OL    H41   49   0.211   0.169   2.235
    35OL    H42   50   0.389   0.170   2.235
    35OL    C05   51   0.300   0.319   2.110
    35OL    H51   52   0.211   0.382   2.109
    35OL    H52   53   0.300   0.255   2.022
    35OL    H53   54   0.389   0.382   2.109
    45OL     HO   55   0.900   1.199   2.801
    45OL     OH   56   0.900   1.142   2.724
    45OL    C01   57   0.900   1.222   2.608
    45OL    H11   58   0.989   1.285   2.608
    45OL    H12   59   0.811   1.285   2.608
    45OL    C02   60   0.900   1.134   2.484
    45OL    H21   61   0.811   1.071   2.484
    45OL    H22   62   0.989   1.071   2.484
    45OL    C03   63   0.900   1.221   2.359
    45OL    H31   64   0.989   1.284   2.359
    45OL    H32   65   0.811   1.284   2.358
    45OL    C04   66   0.900   1.132   2.235
    45OL    H41   67   0.811   1.069   2.235
    45OL    H42   68   0.989   1.070   2.235
    45OL    C05   69   0.900   1.219   2.110
    45OL    H51   70   0.811   1.282   2.109
    45OL    H52   71   0.900   1.155   2.022
    45OL    H53   72   0.989   1.282   2.109
   1.34724   1.55566  10.00000
;
;       File 'mdout.mdp' was generated
;       By user: spoel (291)
;       On host: chagall
;       At date: Mon Dec 15 13:13:06 2003
;

; VARIOUS PREPROCESSING OPTIONS
title                    = test for pentanol
cpp                      = /lib/cpp
include                  =
define                   =

; RUN CONTROL PARAMETERS
integrator               = md
; Start time and timestep in ps
tinit                    = 0
dt                       = 0.001
nsteps                   = 200000
; For exact run continuation or redoing part of a run
init_step                = 0
; mode for center of mass motion removal
comm-mode                = Linear
; number of steps for center of mass motion removal
nstcomm                  = 50

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout                  = 100
nstvout                  = 0
nstfout                  = 0
; Checkpointing helps you continue after crashes
nstcheckpoint            = 0
; Output frequency for energies to log file and energy file
nstlog                   = 10
nstenergy                = 0
; Output frequency and precision for xtc file
nstxtcout                = 0
xtc-precision            = 1000
; This selects the subset of atoms for the xtc file. You can
; select multiple groups. By default all atoms will be written.
xtc-grps                 = 
; Selection of energy groups
energygrps               = 

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist                  = 5
; ns algorithm (simple or grid)
ns_type                  = grid
; Periodic boundary conditions: xyz (default), no (vacuum)
; or full (infinite systems only)
pbc                      = xyz
; nblist cut-off        
rlist                    = 0.55
domain-decomposition     = no

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype              = pme
rcoulomb-switch          = 0
rcoulomb                 = 0.55
; Dielectric constant (DC) for cut-off or DC of reaction field
epsilon-r                = 1
; Method for doing Van der Waals
vdw-type                 = Cut-off
; cut-off lengths       
rvdw-switch              = 0
rvdw                     = 0.55
; Apply long range dispersion corrections for Energy and Pressure
DispCorr                 = EnerPres
; Extension of the potential lookup tables beyond the cut-off
table-extension          = 1
; Spacing for the PME/PPPM FFT grid
fourierspacing           = 0.12
; FFT grid size, when a value is 0 fourierspacing will be used
fourier_nx               = 0
fourier_ny               = 0
fourier_nz               = 0
; EWALD/PME/PPPM parameters
pme_order                = 4
ewald_rtol               = 1e-05
ewald_geometry           = 3dc
epsilon_surface          = 0
optimize_fft             = no

; OPTIONS FOR WEAK COUPLING ALGORITHMS
; Temperature coupling  
Tcoupl                   = nose-hoover
; Groups to couple separately
tc-grps                  = System
; Time constant (ps) and reference temperature (K)
tau_t                    = 0.1
ref_t                    = 10
; Pressure coupling     
Pcoupl                   = no

; SIMULATED ANNEALING  
; Type of annealing for each temperature group (no/single/periodic)
annealing                = single
; Number of time points to use for specifying annealing in each group
annealing_npoints        = 3
; List of times at the annealing points for each group
annealing_time           =  0 10 20
; Temp. at each annealing point, for each group.
annealing_temp           = 10 100 250 

; GENERATE VELOCITIES FOR STARTUP RUN
gen_vel                  = yes
gen_temp                 = 10
gen_seed                 = 1993

; OPTIONS FOR BONDS    
constraints              = none
; Type of constraint algorithm
constraint-algorithm     = Lincs
; Do not constrain the start configuration
unconstrained-start      = no
; Use successive overrelaxation to reduce the number of shake iterations
Shake-SOR                = no
; Relative tolerance of shake
shake-tol                = 1e-04
; Highest order in the expansion of the constraint coupling matrix
lincs-order              = 4
; Number of iterations in the final step of LINCS. 1 is fine for
; normal simulations, but use 2 to conserve energy in NVE runs.
; For energy minimization with constraints it should be 4 to 8.
lincs-iter               = 1
; Lincs will write a warning to the stderr if in one step a bond
; rotates over more degrees than
lincs-warnangle          = 30
; Convert harmonic bonds to morse potentials
morse                    = no


_______________________________________________
gmx-users mailing list    gmx-users@gromacs.org
http://www.gromacs.org/mailman/listinfo/gmx-users
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php

Reply via email to