[gmx-users] Do hydrogen bonds innately exist using mdrun?
Peter Kroon
p.c.kroon at rug.nl
Mon Feb 11 11:21:12 CET 2019
Hi James,
hydrogen bonds are an intrinsically quantum phenomenon, and are more
like a partial covalent bond than anything else. (According to NMR
experiments. Polarisation is transfered across H-bonds).
The good news is that for most purposes MD water "just works". Even
though the exact potential is not present in the model, this gets
compensated by e.g. increased LJ or coulombic interactions. And more
importantly, the rest of the FF (e.g. protein model) is parametrised to
deal with it and compensate further so that in the end the ensemble
sampled during the simulation is physically relevant. This also explains
why some force fields reproduce some behaviour/properties better than
others.
So it really depends on what you need/want.
Peter
On 10-02-19 19:48, James wrote:
> Hi Justin,
>
> Thank you for the information. This is a bit disheartening:
>
> "There is no water model that achieves correct agreement with all
> structural, dynamic, and thermodynamic properties of water"
>
>
> But, I guess that's what you get with force fields rather than using QM
> methods. If you have any familiarity with them, do you think any of the
> reactive force fields are better? (Though even if they are I gather they
> are about 50X slower, which is essentially intractable for many models)
>
> Sincerely,
> James Ryley, PhD, Patent Agent
>
> This communication is confidential and may be subject to legal privilege.
> Nothing contained herein should be construed as legal or patenting advice.
>
>
> On Sat, Feb 9, 2019 at 1:48 AM <
> gromacs.org_gmx-users-request at maillist.sys.kth.se> wrote:
>
>> Send gromacs.org_gmx-users mailing list submissions to
>> gromacs.org_gmx-users at maillist.sys.kth.se
>>
>> To subscribe or unsubscribe via the World Wide Web, visit
>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users
>> or, via email, send a message with subject or body 'help' to
>> gromacs.org_gmx-users-request at maillist.sys.kth.se
>>
>> You can reach the person managing the list at
>> gromacs.org_gmx-users-owner at maillist.sys.kth.se
>>
>> When replying, please edit your Subject line so it is more specific
>> than "Re: Contents of gromacs.org_gmx-users digest..."
>>
>>
>> Today's Topics:
>>
>> 1. Do hydrogen bonds innately exist using mdrun? (James)
>> 2. Re: Do hydrogen bonds innately exist using mdrun? (Justin Lemkul)
>> 3. non-bond (Coulomb, van der Waals) interactions scaled down to
>> 40% of their original value (Jianna Blocchi)
>> 4. Re: Fatal Lambda error (Magnus Lundborg)
>> 5. Simulated Annealing (Shan Jayasinghe)
>>
>>
>> ----------------------------------------------------------------------
>>
>> Message: 1
>> Date: Fri, 8 Feb 2019 14:20:27 -0800
>> From: James <james at ryley.com>
>> To: GROMACS Listserv <gromacs.org_gmx-users at maillist.sys.kth.se>
>> Subject: [gmx-users] Do hydrogen bonds innately exist using mdrun?
>> Message-ID:
>> <
>> CAEhCax-tzvT3n8HaTAaGgVzNAFEsD46L+BscrsZeUWgbHgvMNg at mail.gmail.com>
>> Content-Type: text/plain; charset="UTF-8"
>>
>> Hi,
>>
>> I know that gromacs has an hbond command which will analyze hbonds after a
>> run. But, during a run, without doing something like specifying partial
>> charges, do hbonds exist?
>>
>> It seems to me that they do not, because if I put two water molecules next
>> to each other, they do not assume the configuration they would if hbonds
>> were present, nor do the energies come out as expected (the affinity of the
>> two water molecules is 10X too low, which makes me think that gromacs is
>> using VDW, but ignoring hbonds).
>>
>> I'd love to be wrong and find out that hbonds "just work"! But, if they do
>> not "just work," how do you specify them? I could see mimicking them by
>> tuning partial charges on the donor and acceptor atoms to give the correct
>> affinity, but it seems like that still would not be accurate because hbonds
>> are angle dependent, whereas long range non-bonded forces are not.
>>
>> Sincerely,
>> James Ryley, PhD, Patent Agent
>>
>> This communication is confidential and may be subject to legal privilege.
>> Nothing contained herein should be construed as legal or patenting advice.
>>
>>
>> ------------------------------
>>
>> Message: 2
>> Date: Fri, 8 Feb 2019 19:40:05 -0500
>> From: Justin Lemkul <jalemkul at vt.edu>
>> To: gmx-users at gromacs.org
>> Subject: Re: [gmx-users] Do hydrogen bonds innately exist using mdrun?
>> Message-ID: <cdaa380d-1c06-7b68-e58d-53c34e5a145c at vt.edu>
>> Content-Type: text/plain; charset=utf-8; format=flowed
>>
>>
>>
>> On 2/8/19 5:20 PM, James wrote:
>>> Hi,
>>>
>>> I know that gromacs has an hbond command which will analyze hbonds after
>> a
>>> run. But, during a run, without doing something like specifying partial
>>> charges, do hbonds exist?
>>>
>>> It seems to me that they do not, because if I put two water molecules
>> next
>>> to each other, they do not assume the configuration they would if hbonds
>>> were present, nor do the energies come out as expected (the affinity of
>> the
>>> two water molecules is 10X too low, which makes me think that gromacs is
>>> using VDW, but ignoring hbonds).
>>>
>>> I'd love to be wrong and find out that hbonds "just work"! But, if they
>> do
>>> not "just work," how do you specify them? I could see mimicking them by
>>> tuning partial charges on the donor and acceptor atoms to give the
>> correct
>>> affinity, but it seems like that still would not be accurate because
>> hbonds
>>> are angle dependent, whereas long range non-bonded forces are not.
>> There is no explicit energy function associated with hydrogen bonds.
>> Many force fields used to include such a term, but it was found to be
>> unnecessary with better parametrization. Hydrogen bonds arise due to
>> favorable electrostatic interactions and are entirely a function of the
>> water topology. Don't mess with topologies, because as you "fix" one
>> thing, you break a whole host of others. There is no water model that
>> achieves correct agreement with all structural, dynamic, and
>> thermodynamic properties of water. Most do a good job of hydrogen
>> bonding geometries, RDF, etc.
>>
>> -Justin
>>
>> --
>> ==================================================
>>
>> Justin A. Lemkul, Ph.D.
>> Assistant Professor
>> Office: 301 Fralin Hall
>> Lab: 303 Engel Hall
>>
>> Virginia Tech Department of Biochemistry
>> 340 West Campus Dr.
>> Blacksburg, VA 24061
>>
>> jalemkul at vt.edu | (540) 231-3129
>> http://www.thelemkullab.com
>>
>> ==================================================
>>
>>
>>
>> ------------------------------
>>
>> Message: 3
>> Date: Sat, 9 Feb 2019 09:14:36 +0000 (UTC)
>> From: Jianna Blocchi <jianna_blocchi at yahoo.com>
>> To: "gmx-users at gromacs.org" <gmx-users at gromacs.org>
>> Subject: [gmx-users] non-bond (Coulomb, van der Waals) interactions
>> scaled down to 40% of their original value
>> Message-ID: <1484443140.700236.1549703676420 at mail.yahoo.com>
>> Content-Type: text/plain; charset=UTF-8
>>
>> HiHow can I scale down the?non-bond (Coulomb, van der Waals) interactions
>> to 40% of their original value in OPLS-AA forcefield?Thanks
>>
>> ------------------------------
>>
>> Message: 4
>> Date: Sat, 9 Feb 2019 10:25:29 +0100
>> From: Magnus Lundborg <magnus.lundborg at scilifelab.se>
>> To: gmx-users at gromacs.org
>> Subject: Re: [gmx-users] Fatal Lambda error
>> Message-ID:
>> <CAE7JUpYR677KP=
>> L0CEA62KodeNQnXmajhyFAsEiUAPm+4WT4RA at mail.gmail.com>
>> Content-Type: text/plain; charset="UTF-8"
>>
>> The lambda states are numbered from 0. You only have 21 different states in
>> your mdp file. Therefore your init-lambda-state can range from 0 to 20.
>>
>> Cheers,
>>
>> Magnus
>>
>> Den fre 8 feb. 2019 23:08 skrev Angelina Malagodi <
>> angelinamalagodi at gmail.com>:
>>
>>> Sorry, if the attachments don't load, here are the uploaded texts.
>>> Thank you once again,
>>>
>>> Angelina
>>>
>>>
>>> *Job submission: *
>>>
>>> #!/bin/bash
>>>
>>> # Set some environment variables
>>> FREE_ENERGY=`pwd`
>>> echo "Free energy home directory set to $FREE_ENERGY"
>>> MDP=$FREE_ENERGY/MDP
>>> echo ".mdp files are stored in $MDP"
>>>
>>> # Change to the location of your GROMACS-2018 installation
>>> GMX=/usr/local/gromacs/bin
>>>
>>> for (( i=0; i<22; i++ ))
>>> do
>>> LAMBDA=$i
>>>
>>> # A new directory will be created for each value of lambda and
>>> # at each step in the workflow for maximum organization.
>>>
>>> mkdir Lambda_$LAMBDA
>>> cd Lambda_$LAMBDA
>>>
>>> ##############################
>>> # ENERGY MINIMIZATION STEEP #
>>> ##############################
>>> echo "Starting minimization for lambda = $LAMBDA..."
>>>
>>> mkdir EM
>>> cd EM
>>>
>>> # Iterative calls to grompp and mdrun to run the simulations
>>>
>>> $GMX/gmx grompp -f $MDP/EM/em_steep_$LAMBDA.mdp -c
>>> $FREE_ENERGY/curcumin/solv_curenol_54a7.gro -p
>>> $FREE_ENERGY/curcumin/topol.top -o min$LAMBDA.tpr
>>>
>>> $GMX/gmx mdrun -deffnm min$LAMBDA -ntmpi 1
>>>
>>> sleep 10
>>>
>>>
>>> #####################
>>> # NVT EQUILIBRATION #
>>> #####################
>>> echo "Starting constant volume equilibration..."
>>>
>>> cd ../
>>> mkdir NVT
>>> cd NVT
>>>
>>> $GMX/gmx grompp -f $MDP/NVT/nvt_$LAMBDA.mdp -c ../EM/min$LAMBDA.gro
>> -p
>>> $FREE_ENERGY/curcumin/topol.top -o nvt$LAMBDA.tpr
>>>
>>> $GMX/gmx mdrun -deffnm nvt$LAMBDA -ntmpi 1
>>>
>>> echo "Constant volume equilibration complete."
>>>
>>> sleep 10
>>>
>>> #####################
>>> # NPT EQUILIBRATION #
>>> #####################
>>> echo "Starting constant pressure equilibration..."
>>>
>>> cd ../
>>> mkdir NPT
>>> cd NPT
>>>
>>> $GMX/gmx grompp -f $MDP/NPT/npt_$LAMBDA.mdp -c ../NVT/nvt$LAMBDA.gro
>> -p
>>> $FREE_ENERGY/curcumin/topol.top -t ../NVT/nvt$LAMBDA.cpt -o
>> npt$LAMBDA.tpr
>>> $GMX/gmx mdrun -deffnm npt$LAMBDA -ntmpi 1
>>>
>>> echo "Constant pressure equilibration complete."
>>>
>>> sleep 10
>>>
>>>
>>> #################
>>> # PRODUCTION MD #
>>> #################
>>> echo "Starting production MD simulation..."
>>>
>>> cd ../
>>> mkdir Production_MD
>>> cd Production_MD
>>>
>>> $GMX/gmx grompp -f $MDP/Production_MD/md_$LAMBDA.mdp -c
>>> ../NPT/npt$LAMBDA.gro -p $FREE_ENERGY/curcumin/topol.top -t
>>> ../NPT/npt$LAMBDA.cpt -o md$LAMBDA.tpr
>>>
>>> $GMX/gmx mdrun -deffnm md$LAMBDA -ntmpi 1
>>>
>>> echo "Production MD complete."
>>>
>>> # End
>>> echo "Ending. Job completed for lambda = $LAMBDA"
>>>
>>> cd $FREE_ENERGY
>>> done
>>>
>>> exit;
>>>
>>> *MDP input file:*
>>>
>>> ; Run control
>>> integrator = steep
>>> nsteps = 5000
>>> ; EM criteria and other stuff
>>> emtol = 100
>>> emstep = 0.01
>>> niter = 20
>>> nbfgscorr = 10
>>> ; Output control
>>> nstlog = 1
>>> nstenergy = 1
>>> ; Neighborsearching and short-range nonbonded interactions
>>> cutoff-scheme = verlet
>>> nstlist = 1
>>> ns_type = grid
>>> pbc = xyz
>>> rlist = 1.2
>>> ; Electrostatics
>>> coulombtype = PME
>>> rcoulomb = 1.2
>>> ; van der Waals
>>> vdwtype = cutoff
>>> vdw-modifier = potential-switch
>>> rvdw-switch = 1.0
>>> rvdw = 1.2
>>> ; 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
>>> ; Temperature and pressure coupling are off during EM
>>> tcoupl = no
>>> pcoupl = no
>>> ; Free energy control stuff
>>> free_energy = yes
>>> init_lambda_state = 21
>>> delta_lambda = 0
>>> calc_lambda_neighbors = 1 ; only immediate neighboring windows
>>> couple-moltype = CUR ; 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
>>> ; 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 21
>>> vdw_lambdas = 0.00 0.005 0.01 0.025 0.05 0.075 0.10 0.20
>> 0.30
>>> 0.40 0.50 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.90 0.95 1.00
>>> coul_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
>>> ; 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
>>> ; Not doing simulated temperting here
>>> temperature_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
>>> ; Options for the decoupling
>>> sc-alpha = 0.5
>>> sc-coul = no ; linear interpolation of Coulomb
>> (none
>>> in this case)
>>> sc-power = 1
>>> sc-sigma = 0.3
>>> nstdhdl = 10
>>> ; No velocities during EM
>>> gen_vel = no
>>> ; 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
>>>
>>>
>>>
>>>
>>>
>>> On Fri, Feb 8, 2019 at 2:50 PM Angelina Malagodi <
>>> angelinamalagodi at gmail.com>
>>> wrote:
>>>
>>>> I am trying to calculate the free energy of curcumin in water. My
>> system
>>>> goes from 0- 21 lambdas, but it is not recognizing my final lambda. I
>> am
>>>> sure that I have 21 transitions for vdw_lambdas, coul_lambdas,
>>>> bounded_lambdas, resistant_lambdas, mass_lambdas, and
>>> temperature_lambdas.
>>>> I do not understand why the 21st lambda is not being recognized. My job
>>>> script, input file, and error are attached. Thank you so much for your
>>>> time and help!
>>>>
>>>> Sincerely,
>>>> Angelina Malagodi
>>>>
>>>>
>>>> ERROR 1 [file
>>>> /Users/tyemartin/Documents/Free_Energy/MDP/EM/em_steep_21.mdp]:
>>>>
>>>> initial thermodynamic state 21 does not exist, only goes to 20
>>>>
>>>>
>>>>
>>> --
>>> Gromacs Users mailing list
>>>
>>> * Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>>> posting!
>>>
>>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>
>>> * For (un)subscribe requests visit
>>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>>> send a mail to gmx-users-request at gromacs.org.
>>>
>>
>> ------------------------------
>>
>> Message: 5
>> Date: Sat, 9 Feb 2019 20:47:36 +1100
>> From: Shan Jayasinghe <shanjayasinghe2011 at gmail.com>
>> To: gmx-users at gromacs.org
>> Subject: [gmx-users] Simulated Annealing
>> Message-ID:
>> <
>> CAChmq0h6Zwvgvq3O4TKeuJaSAWS2qEBOGqeLACRHRJro5Wc3Lw at mail.gmail.com>
>> Content-Type: text/plain; charset="UTF-8"
>>
>> Dear Gromacs Users,
>>
>> How do we extend a simulated annealing?
>>
>> Thank you.
>>
>> Best Regards
>> Shan Jayasinghe
>>
>>
>> ------------------------------
>>
>> --
>> Gromacs Users mailing list
>>
>> * Please search the archive at
>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>> posting!
>>
>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>> * For (un)subscribe requests visit
>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>> send a mail to gmx-users-request at gromacs.org.
>>
>> End of gromacs.org_gmx-users Digest, Vol 178, Issue 37
>> ******************************************************
>>
More information about the gromacs.org_gmx-users
mailing list