Energy Analysis¶

pytraj.energy_analysis.esander(traj=None, prmtop=None, igb=8, mm_options=None, qm_options=None, mode=None, dtype='dict', frame_indices=None, top=None)¶

energy decomposition by calling libsander

Parameters:

traj : Trajectory-like or iterables that produce Frame if traj does not hold Topology information, top must be provided prmtop : str or Structure from ParmEd, default=None, optional To avoid any unexpected error, you should always provide original topology filename. If prmtop is None, pytraj will load Topology from traj.top.filename. why do you need to load additional topology filename? Because cpptraj and sander use different Topology object, can not convert from one to another. igb : GB model, default=8 (GB-Neck2) If specify mm_options, this igb input will be ignored mm_options : InputOptions from sander, default=None, optional if mm_options is None, use gas_input with given igb. If mm_options is not None, use this qm_options : InputOptions from sander for QMMM, optional mode : str, default=None, optional if mode=’minimal’, get only ‘bond’, ‘angle’, ‘dihedral’ and ‘total’ energies top : pytraj.Topology or str, default=None, optional only need to specify this top if traj does not hold Topology dtype : str, {‘dict’, ‘dataset’, ‘ndarray’, ‘dataframe’}, default=’dict’ return data type frame_indices : None or 1D array-like, default None if not None, only perform calculation for given frames

Returns:

Dict of energies (to be used with DataFrame) or DatasetList

Notes

This method does not work with pytraj.pmap when you specify mm_options and qm_options. Use pytraj.pmap_mpi with MPI instead.

Work with pytraj.pmap:

pt.pmap(pt.esander, traj, igb=8, dtype='dict')

Will NOT work with pytraj.pmap:

import sander
inp = sander.gas_input(8)
pt.pmap(pt.esander, traj, mm_options=inp, dtype='dict')

Why? Because Python need to pickle each object to send to different cores and Python does not know how to pickle mm_options from sander.gas_input(8).

This works with pytraj.pmap_mpi because pytraj explicitly create mm_options in each core without pickling.

Examples

Examples are adapted from $AMBERHOME/test/sanderapi

>>> import pytraj as pt
>>> # GB energy
>>> traj = pt.datafiles.load_ala3()
>>> traj.n_frames
1
>>> data = pt.esander(traj, igb=8)
>>> data['gb']
array([-92.88577683])
>>> data['bond']
array([ 5.59350521])

>>> # PME
>>> import os
>>> from pytraj.testing import amberhome
>>> import sander
>>> topfile = os.path.join(amberhome, "test/4096wat/prmtop")
>>> rstfile = os.path.join(amberhome, "test/4096wat/eq1.x")
>>> traj = pt.iterload(rstfile, topfile)
>>> options = sander.pme_input()
>>> options.cut = 8.0
>>> edict = pt.esander(traj=traj, mm_options=options)
>>> edict['vdw']
array([ 6028.95167558])

>>> # GB + QMMM
>>> topfile = os.path.join(amberhome, "test/qmmm2/lysine_PM3_qmgb2/prmtop")
>>> rstfile = os.path.join(amberhome, "test/qmmm2/lysine_PM3_qmgb2/lysine.crd")
>>> traj = pt.iterload(rstfile, topfile)

>>> options = sander.gas_input(8)
>>> options.cut = 99.0
>>> options.ifqnt = 1
>>> qm_options = sander.qm_input()
>>> qm_options.iqmatoms[:3] = [8, 9, 10]
>>> qm_options.qm_theory = "PM3"
>>> qm_options.qmcharge = 0
>>> qm_options.qmgb = 2
>>> qm_options.adjust_q = 0

>>> edict = pt.esander(traj=traj, mm_options=options, qm_options=qm_options)
>>> edict['bond']
array([ 0.00160733])
>>> edict['scf']
array([-11.92177575])

>>> # passing options to `pytraj.pmap`: need to pass string
>>> from pytraj.testing import get_fn
>>> fn, tn = get_fn('tz2')
>>> traj = pt.iterload(fn, tn)
>>> inp_str = 'mm_options = sander.pme_input()'
>>> edict = pt.pmap(pt.esander, traj, mm_options=inp_str, n_cores=2)
>>> edict['dihedral']
array([ 126.39307126,  127.0460586 ,  137.26793522,  125.30521069,
        125.25110884,  137.69287326,  125.78280543,  125.14530517,
        118.41540102,  128.73535036])