3

我在 GROMOS54a7 中运行一个简单的苯模拟。我想使用 MDAnalysis 1.0.0 计算每个苯分子质心的 RDF。

这可能吗？我使用 Jupyter Notebook 中的以下代码为 C 分子 g_cc(r) 创建了 rdf：

import MDAnalysis
import numpy as np
%matplotlib inline
import MDAnalysis.analysis.rdf as mda
import matplotlib.pyplot as plt

u = MDAnalysis.Universe("739-c6h6-MolDynamics.tpr","739-c6h6-MolDynamics_good-pbc.xtc")
s1 = u.select_atoms("resid 0 and type CAro")
s2 = u.select_atoms("not (resid 0) and type CAro")
rdf = mda.InterRDF(s1, s2)
rdf.run()

我想获取每个苯分子（每个苯分子在我的模拟中都是一个残基），计算它的 COM 并在上面运行一个类似上面的脚本。有可能做这样的事情吗？

关于 RDF 的一般性问题：我上面使用的方法是否使用我轨迹的每一帧来构造 RDF？我不知道这是否在文档中明确，所以如果这是一个明显的问题，我深表歉意。

感谢您的任何建议！

4

1 回答 1

3

将 CG 基团用作原生原子以重用 MDAnalysis 中的分析工具将很有用。

这是一个模仿 MDAnalysis 组并呈现新positions属性的快速修复。新positions提供几何中心而不是实际位置。我还覆盖了len以传达只有一个珠子用于 CG 元素。

import MDAnalysis as mda
import numpy as np
import MDAnalysis.analysis.rdf
import matplotlib.pyplot as plt

u = mda.Universe('sys_solv.pdb','prod.dcd')

class CG:
    def __init__(self, ag):
        self.ag = ag
        self.universe = self.ag.universe
        self.trajectory = self.ag.universe.trajectory

    @property
    def positions(self):
        return np.array([self.ag.center_of_geometry()])

    def __len__(self):
        return 1

cg_selection = u.select_atoms('resid 1')
cg_atom = CG(cg_selection.atoms)
waters = u.select_atoms('name O')

rdf = MDAnalysis.analysis.rdf.InterRDF(cg_atom, waters)
rdf.run() 
plt.plot(rdf.bins, rdf.rdf)

验证：我为 CG 珠子选择了一个原子，它再现了原始 RDF。

MDAnalysis 使用整个轨迹。您可能会在文档中找到 .run() 函数的启动/停止/步骤参数，该函数允许您缩小要专门使用的帧的范围。

于 2021-02-24T10:40:14.340 回答