Subject: Thermodynamic Integration of Loosening Restraints
Hello,
I would like to determine the free energy change associated with loosening restraints on a small molecule, going from 100 kcal/mol*A^2 to 0. Because this alchemical transformation does not involve a change in the molecule topology, I do not see a straightforward way of doing it in AMBER/pmemd. My best guess at how to accomplish the change was to linearly decrease the restraint_wt parameter for each of my 15 lambda values between 0.0 and 1.0. See below for an example for the first lambda value:
restraint_wt = 100.00, ! Harmonic restraint in kcal/mol*A^2
restraintmask='!:WAT', ! All heavy restrained
icfe = 1, ! Free energy calculations
clambda = 0.0, ! Lambda value for this run
logdvdl = 1, ! Print out dV/dl at end of run
timask1 = '!:WAT', ! Unique atoms at start-state
timask2= '!:WAT', ! Unique atoms at end-state
For my next .in file, restraint_wt would be 95.0 and clambda would be 0.05, and so on. At first I left timask1 and timask2 blank, and I got the error message "ERROR: timask1/2 must match at least one atom for non-softcore run", so I set them both to '*', which threw an error stating that the max number of TI atoms was 500. So, I settled on telling it the solute atoms (!:WAT) were unique to both states. When I run this, the average energies reported at the end of the .out file have "**************" as the value for VDWAALS, EELEC, and VOLUME. Further, VIRIAL is reported as -12582912.2806.
What am I doing wrong? Does AMBER have the capability to use TI to calculate the free energy of going from a restrained small molecule to a flexible one? If not, could I somehow do the calculations "by hand" using data from AMBER simulations?
Thank you,
Helmut Carter
Graduate Student, CUNY
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Feb 22 2024 - 14:00:02 PST
