[gmx-users] Gromacs 4.6 & 4.5.3 qualitative differences & 4.6 instability in polarizable force field vacuum/liquid mixture interface simulations

David van der Spoel spoel at xray.bmc.uu.se
Sat Nov 2 19:52:03 CET 2013


On 2013-11-02 18:38, ploetz wrote:
> Dear Gromacs Users,

Please start a redmine.gromacs.org issue and assign it to me, but try to 
simplify the system as much as possible. You can cut and paste all the 
information to the redmine issue.
>
> I am trying to simulate a system consisting of a vacuum/condensed phase
> interface in which a 6x6x12nm condensed phase region is flanked on both ends
> (in the z-dimension) by a 6x6x12nm vacuum region to form overall box
> dimensions of 6x6x36 nm. The system is a binary liquid mixture of methanol
> (0.125 mole fraction methanol) in water using a polarizable (charge on a
> spring) force field (COS/M methanol and COS/G2 water) at 300K and 1bar. The
> system is stable in Gromacs 4.5.3; however, mdrun gives a segmentation fault
> in Gromacs 4.6 when attempting to do dynamics (energy minimization completes
> with no apparent problems). If I remove the vacuum region, mdrun works. If I
> incrementally add 2 Angstroms to the z-dimension until I reached a vacuum
> region of 34 Angstroms total (17 Angstroms on both sides of the condensed
> phase region) and try to simulate these systems, mdrun works every time.
> When I reach 36 Angstroms, the segmentation fault re-appears. Although not
> the system I am actually interested in, I did some simulations using Gromacs
> 4.6 with the 34 Angstrom vacuum region system and observed an undulating and
> very turbulant "vacuum"/condensed phase interface in which a column of
> water/methanol mixture came out of the condensed phase region to connect the
> two interfaces. Also, the center of mass motion of the system appeared to
> not have been removed. In contrast, using Gromacs 4.5.3, the interface is
> not undulating, but is "calm" and qualitatively planar, no column forms to
> connect the interfaces, and there is no problem with the center of mass
> motion removal. Some of the molecules do enter the "vacuum" region when
> running with 4.5.3, but this appears to be due to the movement of individual
> molecules, not a collective motion of many molecules. This system runs fine
> with a non-polarizable force field in Gromacs 4.6.
>
> I have also compared several properties (using g_energy) of the bulk system
> (no vacuum region) using Gromacs 4.6 and 4.5.3 and they are not the same for
> the polarizable force field, but they are the same for the non-polarizable
> force field. Specifically, with the polarizable force field, the LJ(SR)
> energy is more positive with 4.6, the LJ(LR) energy is more negative with
> 4.6, the Coulomb(SR) energy is more negative with 4.6, the Could. recip.
> energy is more negative with 4.6, the polarization energy is more positive
> with 4.6, the potential energy is more negative in 4.6, the average kinetic
> energy (and temperature) is the same but the fluctuations are greater in
> 4.6, the total energy is more negative in 4.6, the pressure looks fine, and
> the volume looks fine. Where I've noted differences, these are all
> statistically significant differences.
>
> I would like to know if I can just use 4.5.3 and assume the differences
> between the results of 4.5.3 and 4.6 are due to some problem in 4.6. I am
> using all the same input files and commands with both versions, only
> different executables. I ran the regression tests for 4.6 when I installed
> it, and passed them all.
>
> Sincerely,
>
> Elizabeth
>
> -----
> ADDITIONAL DETAILS:
> Below is where the problem appears for the interface system when z-dimension
> of the vacuum region is >= 36 Angstroms total. eq4 is my first attempt at
> dynamics, after three successful energy minimizations (1st: charges screened
> and no bond constraints, 2nd: charges felt but no bond constraints, 3rd:
> charges felt and bond constraints on)
> [ploetz at cluster AddRemainingVacuumBack]$ grompp -f eq4.mdp -c em3.pdb -o
> eq4.tpr -n index.ndx -p sys.top -nice 0
>                           :-)  G  R  O  M  A  C  S  (-:
>
>                    Green Red Orange Magenta Azure Cyan Skyblue
>
>                               :-)  VERSION 4.6  (-:
>
>          Contributions from Mark Abraham, Emile Apol, Rossen Apostolov,
>             Herman J.C. Berendsen, Aldert van Buuren, Pär Bjelkmar,
>       Rudi van Drunen, Anton Feenstra, Gerrit Groenhof, Christoph Junghans,
>          Peter Kasson, Carsten Kutzner, Per Larsson, Pieter Meulenhoff,
>             Teemu Murtola, Szilard Pall, Sander Pronk, Roland Schulz,
>                  Michael Shirts, Alfons Sijbers, Peter Tieleman,
>
>                 Berk Hess, David van der Spoel, and Erik Lindahl.
>
>         Copyright (c) 1991-2000, University of Groningen, The Netherlands.
>           Copyright (c) 2001-2012,2013, The GROMACS development team at
>          Uppsala University & The Royal Institute of Technology, Sweden.
>              check out http://www.gromacs.org for more information.
>
>           This program is free software; you can redistribute it and/or
>         modify it under the terms of the GNU Lesser General Public License
>          as published by the Free Software Foundation; either version 2.1
>               of the License, or (at your option) any later version.
>
>                                  :-)  grompp  (-:
>
> Option     Filename  Type         Description
> ------------------------------------------------------------
>    -f        eq4.mdp  Input        grompp input file with MD parameters
>   -po      mdout.mdp  Output       grompp input file with MD parameters
>    -c        em3.pdb  Input        Structure file: gro g96 pdb tpr etc.
>    -r       conf.gro  Input, Opt.  Structure file: gro g96 pdb tpr etc.
>   -rb       conf.gro  Input, Opt.  Structure file: gro g96 pdb tpr etc.
>    -n    index.ndx  Input, Opt!  Index file
>    -p     sys.top  Input        Topology file
>   -pp  processed.top  Output, Opt. Topology file
>    -o        eq4.tpr  Output       Run input file: tpr tpb tpa
>    -t       traj.trr  Input, Opt.  Full precision trajectory: trr trj cpt
>    -e       ener.edr  Input, Opt.  Energy file
> -ref     rotref.trr  In/Out, Opt. Full precision trajectory: trr trj cpt
>
> Option       Type   Value   Description
> ------------------------------------------------------
> -[no]h       bool   no      Print help info and quit
> -[no]version bool   no      Print version info and quit
> -nice        int    0       Set the nicelevel
> -[no]v       bool   no      Be loud and noisy
> -time        real   -1      Take frame at or first after this time.
> -[no]rmvsbds bool   yes     Remove constant bonded interactions with virtual
>                              sites
> -maxwarn     int    0       Number of allowed warnings during input
>                              processing. Not for normal use and may generate
>                              unstable systems
> -[no]zero    bool   no      Set parameters for bonded interactions without
>                              defaults to zero instead of generating an error
> -[no]renum   bool   yes     Renumber atomtypes and minimize number of
>                              atomtypes
>
> Ignoring obsolete mdp entry 'cpp'
> Ignoring obsolete mdp entry 'domain-decomposition'
> Ignoring obsolete mdp entry 'andersen_seed'
> Ignoring obsolete mdp entry 'nstcheckpoint'
> Replacing old mdp entry 'unconstrained-start' by 'continuation'
>
> NOTE 1 [file eq4.mdp]:
>    The Berendsen thermostat does not generate the correct kinetic energy
>    distribution. You might want to consider using the V-rescale thermostat.
>
> Generated 80 of the 91 non-bonded parameter combinations
> Excluding 2 bonded neighbours molecule type 'COSM'
> turning all bonds into constraints...
> Excluding 2 bonded neighbours molecule type 'COSG2'
> turning all bonds into constraints...
> Velocities were taken from a Maxwell distribution at 300.15 K
> Cleaning up constraints and constant bonded interactions with virtual sites
> Number of degrees of freedom in T-Coupling group MOH is 8243.62
> Number of degrees of freedom in T-Coupling group SOL is 57717.38
> Largest charge group radii for Van der Waals: 0.118, 0.118 nm
> Largest charge group radii for Coulomb:       0.118, 0.118 nm
> Calculating fourier grid dimensions for X Y Z
> Using a fourier grid of 48x48x288, spacing 0.120 0.120 0.120
> Estimate for the relative computational load of the PME mesh part: 0.37
> This run will generate roughly 622 Mb of data
>
> There was 1 note
>
> gcq#228: "Bailed Out Of Edge Synchronization After 10,000 Iterations"
> (X/Motif)
>
> [ploetz at cluster AddRemainingVacuumBack]$ more eq4.mdp
> ; VARIOUS PREPROCESSING OPTIONS
> title                    =
> ; Preprocessor - specify a full path if necessary.
> cpp                      = /lib/cpp
> include                  =
> define                   =
>
> ; RUN CONTROL PARAMETERS
> integrator               = md
> ; Start time and timestep in ps
> tinit                    = 0
> dt                       = 0.002
> nsteps                   = 500000
> ; 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                  = 500
> ; group(s) for center of mass motion removal
> comm-grps                =
>
> ; OUTPUT CONTROL OPTIONS
> ; Output frequency for coords (x), velocities (v) and forces (f)
> nstxout                  = 500
> nstvout                  = 0
> nstfout                  = 0
> ; Checkpointing helps you continue after crashes
> nstcheckpoint            = 5000
> ; Output frequency for energies to log file and energy file
> nstlog                   = 500
> nstenergy                = 500
> ; 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               = MOH SOL
>
> ; NEIGHBORSEARCHING PARAMETERS
> ; nblist update frequency
> nstlist                  = 10
> ; 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                    = 1.0
> domain-decomposition     = no
>
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> ; Method for doing electrostatics
> coulombtype              = PME
> rcoulomb-switch          = 0
> rcoulomb                 = 1.0
> ; Relative dielectric constant for the medium and the reaction field
> epsilon-r                = 1
> epsilon_rf               = 1
> ; Method for doing Van der Waals
> vdw-type                 = Cut-off
> ; cut-off lengths
> rvdw-switch              = 0
> rvdw                     = 1.5
> ; Apply long range dispersion corrections for Energy and Pressure
> DispCorr                 = No
> ; Extension of the potential lookup tables beyond the cut-off
> table-extension          = 1
> ; Seperate tables between energy group pairs
> energygrp_table          =
> ; 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           = 3d
> epsilon_surface          = 0
> optimize_fft             = no
>
> ; IMPLICIT SOLVENT (for use with Generalized Born electrostatics)
> implicit_solvent         = No
>
> ; OPTIONS FOR WEAK COUPLING ALGORITHMS
> ; Temperature coupling
> tcoupl                   = berendsen
> ; Groups to couple separately
> tc-grps                  = MOH SOL
> ; Time constant (ps) and reference temperature (K)
> tau-t                    = 0.1 0.1
> ref-t                    = 300.15 300.15
> ; Pressure coupling
> Pcoupl                   = no
> Pcoupltype               = isotropic
> ; Time constant (ps), compressibility (1/bar) and reference P (bar)
> tau-p                    = 0.5
> compressibility          = 4.5e-5
> ref-p                    = 1
> ; Random seed for Andersen thermostat
> andersen_seed            = 815131
>
> ; GENERATE VELOCITIES FOR STARTUP RUN
> gen-vel                  = yes
> gen-temp                 = 300.15
> gen-seed                 = 173529
>
> ; OPTIONS FOR BONDS
> constraints              = all-bonds
> ; 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                = 0.0001
> ; 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               = 4
> ; 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
>
> ; ENERGY GROUP EXCLUSIONS
> ; Pairs of energy groups for which all non-bonded interactions are excluded
> energygrp_excl           =
>
> [ploetz at cluster AddRemainingVacuumBack]$
>
> [ploetz at cluster AddRemainingVacuumBack]$ mdrun -v -s eq4.tpr -c eq4.pdb -g
> eq4.log -o eq4.trr -nice 0
>                           :-)  G  R  O  M  A  C  S  (-:
>
>                 GRoups of Organic Molecules in ACtion for Science
>
>                               :-)  VERSION 4.6  (-:
>
>          Contributions from Mark Abraham, Emile Apol, Rossen Apostolov,
>             Herman J.C. Berendsen, Aldert van Buuren, Pär Bjelkmar,
>       Rudi van Drunen, Anton Feenstra, Gerrit Groenhof, Christoph Junghans,
>          Peter Kasson, Carsten Kutzner, Per Larsson, Pieter Meulenhoff,
>             Teemu Murtola, Szilard Pall, Sander Pronk, Roland Schulz,
>                  Michael Shirts, Alfons Sijbers, Peter Tieleman,
>
>                 Berk Hess, David van der Spoel, and Erik Lindahl.
>
>         Copyright (c) 1991-2000, University of Groningen, The Netherlands.
>           Copyright (c) 2001-2012,2013, The GROMACS development team at
>          Uppsala University & The Royal Institute of Technology, Sweden.
>              check out http://www.gromacs.org for more information.
>
>           This program is free software; you can redistribute it and/or
>         modify it under the terms of the GNU Lesser General Public License
>          as published by the Free Software Foundation; either version 2.1
>               of the License, or (at your option) any later version.
>
>                                  :-)  mdrun  (-:
>
> Option     Filename  Type         Description
> ------------------------------------------------------------
>    -s        eq4.tpr  Input        Run input file: tpr tpb tpa
>    -o        eq4.trr  Output       Full precision trajectory: trr trj cpt
>    -x       traj.xtc  Output, Opt. Compressed trajectory (portable xdr
> format)
> -cpi      state.cpt  Input, Opt.  Checkpoint file
> -cpo      state.cpt  Output, Opt. Checkpoint file
>    -c        eq4.pdb  Output       Structure file: gro g96 pdb etc.
>    -e       ener.edr  Output       Energy file
>    -g        eq4.log  Output       Log file
> -dhdl      dhdl.xvg  Output, Opt. xvgr/xmgr file
> -field    field.xvg  Output, Opt. xvgr/xmgr file
> -table    table.xvg  Input, Opt.  xvgr/xmgr file
> -tabletf    tabletf.xvg  Input, Opt.  xvgr/xmgr file
> -tablep  tablep.xvg  Input, Opt.  xvgr/xmgr file
> -tableb   table.xvg  Input, Opt.  xvgr/xmgr file
> -rerun    rerun.xtc  Input, Opt.  Trajectory: xtc trr trj gro g96 pdb cpt
> -tpi        tpi.xvg  Output, Opt. xvgr/xmgr file
> -tpid   tpidist.xvg  Output, Opt. xvgr/xmgr file
>   -ei        sam.edi  Input, Opt.  ED sampling input
>   -eo      edsam.xvg  Output, Opt. xvgr/xmgr file
>    -j       wham.gct  Input, Opt.  General coupling stuff
>   -jo        bam.gct  Output, Opt. General coupling stuff
> -ffout      gct.xvg  Output, Opt. xvgr/xmgr file
> -devout   deviatie.xvg  Output, Opt. xvgr/xmgr file
> -runav  runaver.xvg  Output, Opt. xvgr/xmgr file
>   -px      pullx.xvg  Output, Opt. xvgr/xmgr file
>   -pf      pullf.xvg  Output, Opt. xvgr/xmgr file
>   -ro   rotation.xvg  Output, Opt. xvgr/xmgr file
>   -ra  rotangles.log  Output, Opt. Log file
>   -rs   rotslabs.log  Output, Opt. Log file
>   -rt  rottorque.log  Output, Opt. Log file
> -mtx         nm.mtx  Output, Opt. Hessian matrix
>   -dn     dipole.ndx  Output, Opt. Index file
> -multidir    rundir  Input, Opt., Mult. Run directory
> -membed  membed.dat  Input, Opt.  Generic data file
>   -mp     membed.top  Input, Opt.  Topology file
>   -mn     membed.ndx  Input, Opt.  Index file
>
> Option       Type   Value   Description
> ------------------------------------------------------
> -[no]h       bool   no      Print help info and quit
> -[no]version bool   no      Print version info and quit
> -nice        int    0       Set the nicelevel
> -deffnm      string         Set the default filename for all file options
> -xvg         enum   xmgrace  xvg plot formatting: xmgrace, xmgr or none
> -[no]pd      bool   no      Use particle decompostion
> -dd          vector 0 0 0   Domain decomposition grid, 0 is optimize
> -ddorder     enum   interleave  DD node order: interleave, pp_pme or
> cartesian
> -npme        int    -1      Number of separate nodes to be used for PME, -1
>                              is guess
> -nt          int    0       Total number of threads to start (0 is guess)
> -ntmpi       int    0       Number of thread-MPI threads to start (0 is
> guess)
> -ntomp       int    0       Number of OpenMP threads per MPI process/thread
>                              to start (0 is guess)
> -ntomp_pme   int    0       Number of OpenMP threads per MPI process/thread
>                              to start (0 is -ntomp)
> -pin         enum   auto    Fix threads (or processes) to specific cores:
>                              auto, on or off
> -pinoffset   int    0       The starting logical core number for pinning to
>                              cores; used to avoid pinning threads from
>                              different mdrun instances to the same core
> -pinstride   int    0       Pinning distance in logical cores for threads,
>                              use 0 to minimize the number of threads per
>                              physical core
> -gpu_id      string         List of GPU id's to use
> -[no]ddcheck bool   yes     Check for all bonded interactions with DD
> -rdd         real   0       The maximum distance for bonded interactions
> with
>                              DD (nm), 0 is determine from initial coordinates
> -rcon        real   0       Maximum distance for P-LINCS (nm), 0 is estimate
> -dlb         enum   auto    Dynamic load balancing (with DD): auto, no or
> yes
> -dds         real   0.8     Minimum allowed dlb scaling of the DD cell size
> -gcom        int    -1      Global communication frequency
> -nb          enum   auto    Calculate non-bonded interactions on: auto, cpu,
>                              gpu or gpu_cpu
> -[no]tunepme bool   yes     Optimize PME load between PP/PME nodes or
> GPU/CPU
> -[no]testverlet bool   no      Test the Verlet non-bonded scheme
> -[no]v       bool   yes     Be loud and noisy
> -[no]compact bool   yes     Write a compact log file
> -[no]seppot  bool   no      Write separate V and dVdl terms for each
>                              interaction type and node to the log file(s)
> -pforce      real   -1      Print all forces larger than this (kJ/mol nm)
> -[no]reprod  bool   no      Try to avoid optimizations that affect binary
>                              reproducibility
> -cpt         real   15      Checkpoint interval (minutes)
> -[no]cpnum   bool   no      Keep and number checkpoint files
> -[no]append  bool   yes     Append to previous output files when continuing
>                              from checkpoint instead of adding the simulation
>                              part number to all file names
> -nsteps      int    -2      Run this number of steps, overrides .mdp file
>                              option
> -maxh        real   -1      Terminate after 0.99 times this time (hours)
> -multi       int    0       Do multiple simulations in parallel
> -replex      int    0       Attempt replica exchange periodically with this
>                              period (steps)
> -nex         int    0       Number of random exchanges to carry out each
>                              exchange interval (N^3 is one suggestion).  -nex
>                              zero or not specified gives neighbor replica
>                              exchange.
> -reseed      int    -1      Seed for replica exchange, -1 is generate a seed
> -[no]ionize  bool   no      Do a simulation including the effect of an X-Ray
>                              bombardment on your system
>
> Reading file eq4.tpr, VERSION 4.6 (single precision)
> Using 12 MPI threads
> starting mdrun 'COSM in water x_COSM 0.125'
> 500000 steps,   1000.0 ps.
> Segmentation fault
> [ploetz at cluster AddRemainingVacuumBack]$
>
> [ploetz at cluster AddRemainingVacuumBack]$ more sys.top
> #include "polff.itp"
> #include "cosm.itp"
> #include "cosg2.itp"
> [ system ]
> ; Name
> COSM in water x_COSM 0.125
> [ molecules ]
> ; Compound        #mols
> COSM                1374
> COSG2                9620
>
>
> [ploetz at cluster AddRemainingVacuumBack]$ more polff.itp
>
> ;Polarizable force field atom types
> ;Includes types for COS/G2 (water) COS/M (MeOH 1-site) and CPC (MeOH 2-site)
> [ defaults ]
> LJ      Geometric
> [ atomtypes ]
> ;name        mass      charge   ptype   c6      c12
>     WO    15.99940       0.0     A       0.0     0.0
>     WH     1.00800       0.0     A       0.0     0.0
>     WD     0.00000       0.0     D       0.0     0.0
>     CSMO  15.99940       0.0     A       0.0     0.0
>     CSMH   1.00800       0.0     A       0.0     0.0
>     CSMC  15.03500       0.0     A       0.0     0.0
>     CPCO  15.99940       0.0     A       0.0     0.0
>     CPCH   1.00800       0.0     A       0.0     0.0
>     CPCC  15.03500       0.0     A       0.0     0.0
>     WS     0.00000       0.0     S       0.0     0.0
>     CSMS   0.00000       0.0     S       0.0     0.0
>     CPCSC  0.00000       0.0     S       0.0     0.0
>     CPCSO  0.00000       0.0     S       0.0     0.0
> [ nonbond_params ]
> ;NB: mix of COSM and CPC not included
> ;NB: includes special COSMC - COS/G2 C12  parameters
> WO      WO      1       3.24434e-03          3.45765e-06
> WO      CSMO    1       2.71125e-03          3.00305e-06
> WO      CSMC    1       5.36612e-03          7.71936e-06
> WO      CPCO    1       2.68847e-03          2.74273e-06
> WO      CPCC    1       6.35721e-03         10.57112e-06
> CSMO    CSMO    1       2.26576e-03          2.60823e-06
> CSMO    CSMC    1       4.48440e-03          7.10600e-06
> CSMC    CSMC    1       8.87552e-03         19.36000e-06
> CPCO    CPCO    1       2.22784e-03          2.17563e-06
> CPCO    CPCC    1       5.26799e-03          8.38538e-06
> CPCC    CPCC    1      12.45679e-03         32.31923e-06
> [ploetz at cluster AddRemainingVacuumBack]$
>
> [ploetz at cluster AddRemainingVacuumBack]$  more cosm.itp
>
> ;  Topology for COS/M methanol model
> ;  Yu et al. J., Comput. Chem. 27 (1494-1504), 2006
> [ moleculetype ]
> ; molname       nrexcl
> COSM            2
> [ atoms ]
> ; id    at type res nr  residu name     at name         cg nr   charge
> 1       CSMO    1       MOH             OM              1       7.470
> 2       CSMH    1       MOH             HM              1       0.360
> 3       CSMC    1       MOH             CM              1       0.170
> 4       CSMS    1       MOH            SOM              1      -8.000
> [ polarization ]
> ;Center Shell  funct       alpha (nm^3)
>   1      4       1       0.001320
> [ constraints ]
> ; ai    aj   funct   length
> 1       2      2     0.1000
> 1       3      2     0.1430
> 2       3      2     0.1988
> [ exclusions ]
> ; iatom excluded from interaction with i
> 1       2       3       4
> 2       1       3       4
> 3       1       2       4
> 4       1       2       3
> #ifdef POSRES
> [ position_restraints ]
> ; iatom type    fx      fy      fz
> 1       1       100     100     100
> 3       1       100     100     100
> #endif
> [ploetz at cluster AddRemainingVacuumBack]$
> [ploetz at cluster AddRemainingVacuumBack]$ more cosg2.itp
>
> ;  Topology for COS/G2 water model
> ;  Yu and Van Gunsteren, J. Chem. Phys. 121 (9549-9564), 2004
> [ moleculetype ]
> ; molname       nrexcl
> COSG2           2
> [ atoms ]
> ; id    at type res nr  residu name     at name         cg nr   charge
> 1       WO      1       SOL             OW              1       0.0000
> 2       WH      1       SOL             HW1             1       0.5265
> 3       WH      1       SOL             HW2             1       0.5265
> 4       WD      1       SOL             MW              1       6.9470
> 5       WS      1       SOL             SW              1      -8.0000
> [ polarization ]
> ;Center Shell   funct       alpha (nm^3)
>   4      5       1       0.001255
> [ settles ]
> ; i     funct   dOH     dHH
> 1       1       0.09572 0.15139
> [ dummies3 ]
> ; The position of the dummies is computed as follows:
> ;
> ;               O
> ;
> ;               D
> ;
> ;       H               H
> ;
> ; 2 * b = distance (OD) / [ cos (angle(DOH))    * distance (OH) ]
> ;         0.022 nm      / [ cos (104.52 / 2 deg) * 0.09572 nm   ]
> ;         0.022 nm      / 0.058588228 nm
> ; Dummy pos x4 = x1 + a*(x2-x1) + b*(x3-X1)
> ; Dummy from                    funct   a               b
> 4       1       2       3       1       0.187751028     0.187751028
> [ exclusions ]
> ; iatom excluded from interaction with i
> 1       2       3       4       5
> 2       1       3       4       5
> 3       1       2       4       5
> 4       1       2       3       5
> 5       1       2       3       4
> #ifdef POSRES
> [ position_restraints ]
> ; iatom type    fx      fy      fz
> 1       1       100     100     100
> #endif
> [ploetz at cluster AddRemainingVacuumBack]$
>
> --
> View this message in context: http://gromacs.5086.x6.nabble.com/Gromacs-4-6-4-5-3-qualitative-differences-4-6-instability-in-polarizable-force-field-vacuum-liquid-ms-tp5012174.html
> Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
>


-- 
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205.
spoel at xray.bmc.uu.se    http://folding.bmc.uu.se



More information about the gromacs.org_gmx-users mailing list