Dear Wes, Thanks for clarifying, I did as you told and it worked. The resulting dG_hyd of NH4+ is -62.5 kJ/mol, while the experimental value is -285 kJ/mol (10.1039/FT9918702995). I am using OPLS and SPC/E water. Energy values of the last lambdas showed higher error bars, so now I am refining the last part by increasing the number of lambda points: fep-lambdas = 0.0 0.2 0.4 0.6 0.8 0.85 0.9 0.92 0.94 0.96 0.98 1.0 I do not expect big differences. I could also try (i) other water model, (ii) bigger box with more NH4+ ions (same concentration) or (ii-b) increase NH4+ concentration, (iii) adding anions (like Cl-). Regarding these options, I am not sure yet if they have sense (meaning that there may be another more suitable option, while some of the possible tests I wrote may not result in significative improvements), so if there are any comments regarding these options they are welcome.
Also, I do not understand why the Nitrogen partial charge is negative, shouldnt be positive? I choose OPLS since it was parameterized for liquid simulations, but I found in literature that usually there are problems for correctly simulating dG_hyd of amines. atom charge mass N00 -0.5183 14.0070 H01 0.3796 1.0080 H02 0.3796 1.0080 H03 0.3796 1.0080 H04 0.3796 1.0080 Maybe the best option is to change FF (does someone know about an acceptable FF -preferably all_atoms- for amines/ ammonium derivatives R4N+) ? Best regards Javier Date: Mon, 19 Mar 2018 13:50:19 -0400 > From: Wes Barnett <w.barn...@columbia.edu> > To: gmx-us...@gromacs.org > Subject: Re: [gmx-users] Free Energy solvation generate multiple > dhdl.xvg output files > Message-ID: > <CALkN6e3i4ity-kZOf5c7eZO=5jumdrBPz620LQK_1nr2YSExxw@mail. > gmail.com> > Content-Type: text/plain; charset="UTF-8" > > It sounds like you are only running one free energy simulation. You need to > run one simulation for each lambda (changing "init-lambda-state" in your > mdp file) and name the dhdl xvg file a different name for each one. > > On Mon, Mar 19, 2018 at 8:30 AM, Javier Luque Di Salvo < > jluquedisa...@gmail.com> wrote: > > > Dear Gromacs users, > > I am trying to run Free Energy calculations of a single ion in water (to > > calculate its solvation energy). Since this is a new topic for me, I am > > following the tutorial 'Solvation free energy of ethanol' (easily find in > > google). The three steps work well without erros (1-energy minimization, > > 2-equilibration, 3-lambda (de)coupling) but I am not obtaining the > separate > > dhdl.xvg files after the step 3, only a single output file run.xvg. The > > reason for this is that in the tutorial, a script named mklambdas.sh is > > used to generate such separate files in separate folders. Since I would > > like to keep track of the .mdp strings and understand how they work, I am > > not using this script, instead I am trying to properly define mdp > commands: > > > > _________________________________________ > > ; Free energy variables > > free-energy = yes > > couple-moltype = NH4 > > couple-lambda0 = none > > couple-lambda1 = vdwq > > couple-intramol = no > > init-lambda = -1 > > init-lambda-state = 0 > > delta-lambda = 0 > > nstdhdl = 50 > > fep-lambdas = 0.0 0.2 0.4 0.6 0.8 0.9 1.0 > > mass-lambdas = > > coul-lambdas = > > vdw-lambdas = > > bonded-lambdas = > > restraint-lambdas = > > temperature-lambdas = > > calc-lambda-neighbors = 1 > > init-lambda-weights = > > dhdl-print-energy = no > > sc-alpha = 1.0 > > sc-power = 1 > > sc-r-power = 6 > > sc-sigma = 0.3 > > sc-coul = no > > separate-dhdl-file = yes > > dhdl-derivatives = yes > > dh_hist_size = 0 > > dh_hist_spacing = 0.1 > > _________________________________________ > > > > When trying to use "gmx bar -f run.xvg -g run.edr -o bar.xvg" I get this > > error: > > _______________________________________________________ > > run.xvg: Ignoring set 'pV (kJ/mol)'. > > run.xvg: 0.0 - 200.0; lambda = 0 > > dH/dl & foreign lambdas: > > dH/dl (fep-lambda) (2001 pts) > > delta H to 0 (2001 pts) > > delta H to 0.2 (2001 pts) > > > > Opened run.edr as single precision energy file > > Reading energy frame 0 time 0.000 > > ------------------------------------------------------- > > Program: gmx bar, version 2016.4 > > Source file: src/gromacs/gmxana/gmx_bar.cpp (line 3286) > > Fatal error: > > Did not find delta H information in file run.edr > > ____________________________________________________ > > > > This is the run.xvg output file: > > > > ____________________________________________________ > > # Command line: > > # gmx mdrun -ntmpi 2 -ntomp 8 -pin on -v -deffnm run > > # gmx mdrun is part of G R O M A C S: > > # > > # S C A M O R G > > # > > @ title "dH/d\xl\f{} and \xD\f{}H" > > @ xaxis label "Time (ps)" > > @ yaxis label "dH/d\xl\f{} and \xD\f{}H (kJ/mol [\xl\f{}]\S-1\N)" > > @TYPE xy > > @ subtitle "T = 300 (K) \xl\f{} state 0: fep-lambda = 0.0000" > > @ view 0.15, 0.15, 0.75, 0.85 > > @ legend on > > @ legend box on > > @ legend loctype view > > @ legend 0.78, 0.8 > > @ legend length 2 > > @ s0 legend "dH/d\xl\f{} fep-lambda = 0.0000" > > @ s1 legend "\xD\f{}H \xl\f{} to 0.0000" > > @ s2 legend "\xD\f{}H \xl\f{} to 0.2000" > > @ s3 legend "pV (kJ/mol)" > > 0.0000 -530.14746 0.0000000 -106.84067 9.7838984 > > 0.1000 -177.17302 0.0000000 -33.945797 9.7844772 > > 0.2000 -158.46501 0.0000000 -30.541293 9.7854595 > > 0.3000 -143.72632 0.0000000 -27.669909 9.7857714 > > (...) > > > > ____________________________________________________ > > > > > > Could you guide me through a proper string definition in the mdp options > to > > get separate files, in order to use the BAR algorithm afterwards? > > > > Best regards > > > > *Javier* > > -- > > 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-requ...@gromacs.org. > > > > > > -- > James "Wes" Barnett > Postdoctoral Research Scientist > Department of Chemical Engineering > Kumar Research Group <http://www.columbia.edu/cu/kumargroup/> > Columbia University > w.barn...@columbia.edu > http://wbarnett.us > > > ------------------------------ > > -- > 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-requ...@gromacs.org. > > End of gromacs.org_gmx-users Digest, Vol 167, Issue 93 > ****************************************************** > -- ________________________________________ *Javier Luque Di Salvo* Dipartamento di Ingegneria Chimica Universitá Degli Studi di Palerm*o* -- 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-requ...@gromacs.org.