Energy Decomposition

pytraj.analysis.energy_analysis.esander(traj=None, prmtop=None, igb=8, mm_options=None, qm_options=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

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])
pytraj.analysis.energy_analysis.lie(traj, mask, options='', dtype='dict', frame_indices=None, top=None)

LIE