Hi Jason,
Thanks for the suggestion. I agree I could keep the water and use the sander API to compute the energies.
- Lee-Ping
> On May 25, 2015, at 12:12 PM, Jason Swails <jason.swails.gmail.com> wrote:
>
> On Mon, May 25, 2015 at 12:44 PM, Lee-Ping Wang <leeping.stanford.edu>
> wrote:
>
>> Hi everyone,
>>
>> I’ve been running some replica exchange simulations in Amber. I’ve got a
>> simple script for cpptraj (adapted from the tutorial at
>> http://ambermd.org/tutorials/advanced/tutorial7/ <
>> http://ambermd.org/tutorials/advanced/tutorial7/> ) for processing the
>> output trajectories, resulting in “temperature-based” trajectories where
>> the temperature is constant but the conformations jump around. The
>> simulations are running very nicely; thanks for writing a great software
>> package.
>>
>> My question is how I could get potential energies to “go with” these
>> trajectories? I would like to use these potential energies as input to a
>> WHAM / MBAR calculation to improve the estimation of thermodynamic
>> properties at each temperature.
>>
>
> The standard way is to feed those trajectories into sander, set imin=5 to
> process an input trajectory, and get energies that way. But I'll provide
> an alternative I think you may find more amenabl
> e:
>
> Use the sander API to compute single-point energies (and assuming you're
> going to turn around and feed this directly into pymbar, you can do
> everything you want in a single script). You can do something like this:
>
> import sander
> import chemistry as chem
>
> parm = chem.amber.AmberParm('prmtop', 'inpcrd')
> traj = chem.load_file('your_trajectory') # mdtraj or pytraj also works to
> read trajectories
>
> input = sander.pme_input() # or use "gas_input(NUM)" where NUM is the igb
> model for implicit solvent
> input.cut = 9.0 # struct members have same name as mdin variables
>
> with sander.setup('prmtop', parm.coords, parm.box, input):
> for i in range(traj.frame):
> sander.set_positions(traj.coordinates(i))
> ene, frc = sander.energy_forces()
> # Go ahead and store your potential energy in a numpy array you set
> up
>
> You can also use OpenMM instead of the sander Python API if you're more
> comfortable with that.
>
>
>
>> My simple script is attached below.
>>
>> Thanks,
>>
>> - Lee-Ping
>>
>> parm ../prmtop
>> ensemble ../remd.00/netcdf.00
>> ensemble ../remd.01/netcdf.00
>> ensemble ../remd.02/netcdf.00
>> strip :WAT
>>
>
> Obviously you can't calculate the full energies if you strip waters, but I
> figured I'd point this out anyway.
>
> All the best,
> Jason
>
> --
> Jason M. Swails
> BioMaPS,
> Rutgers University
> Postdoctoral Researcher
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon May 25 2015 - 14:00:02 PDT
