首页 > 解决方案 > 使用 pyscf 可视化电子密度?

问题描述

我找到了以下代码示例,它使用密度泛函理论来计算电子密度 rho

#!/usr/bin/env python
#
# Author: Qiming Sun <osirpt.sun@gmail.com>
#

import numpy
from pyscf import lib
from pyscf.dft import numint, gen_grid

"""
Gaussian cube file format
"""


def density(mol, outfile, dm, nx=80, ny=80, nz=80):
coord = mol.atom_coords()
box = numpy.max(coord, axis=0) - numpy.min(coord, axis=0) + 4
boxorig = numpy.min(coord, axis=0) - 2
xs = numpy.arange(nx) * (box[0] / nx)
ys = numpy.arange(ny) * (box[1] / ny)
zs = numpy.arange(nz) * (box[2] / nz)
coords = lib.cartesian_prod([xs, ys, zs])
coords = numpy.asarray(coords, order="C") - (-boxorig)

nao = mol.nao_nr()
ngrids = nx * ny * nz
blksize = min(200, ngrids)
rho = numpy.empty(ngrids)
for ip0, ip1 in gen_grid.prange(0, ngrids, blksize):
    ao = numint.eval_ao(mol, coords[ip0:ip1])
    rho[ip0:ip1] = numint.eval_rho(mol, ao, dm)
rho = rho.reshape(nx, ny, nz)

with open(outfile, "w") as f:
    f.write("Density in real space\n")
    f.write("Comment line\n")
    f.write("]" % mol.natm)
    f.write(" .8f .8f .8f\n" % tuple(boxorig.tolist()))
    f.write("] .8f .8f .8f\n" % (nx, xs[1], 0, 0))
    f.write("] .8f .8f .8f\n" % (ny, 0, ys[1], 0))
    f.write("] .8f .8f .8f\n" % (nz, 0, 0, zs[1]))
    for ia in range(mol.natm):
        chg = mol.atom_charge(ia)
        f.write("%5d %f" % (chg, chg))
        f.write(" .8f .8f .8f\n" % tuple(coord[ia]))
    fmt = " .8e" * nz + "\n"
    for ix in range(nx):
        for iy in range(ny):
            f.write(fmt % tuple(rho[ix, iy].tolist()))

if name == " main ": from pyscf import gto, scf from pyscf.tools import cubegen

mol = gto.M(atom="H 0 0 0; H 0 0 1")
mf = scf.RHF(mol)
mf.scf()
cubegen.density(mol, "h2.cube", mf.make_rdm1())

然而,我想知道这个输出 rho 的性质。因此,如果我要可视化这个电子密度,我将如何绘制它?有人熟悉这个包吗?它似乎是一个包含三列 \rho 的数组,每列长度为 80。但是一个点对应于什么?

标签: visualizationphysicsdensity-plotchemistry

解决方案


推荐阅读