Skip to content
python
from importlib.metadata import version
print(version('pywfn'))
from datetime import datetime
print(datetime.now())
1.0.18
2026-07-03 23:14:20.407911

空间性质

计算分子空间中某些点的性质,全都是空间函数v=f(x,y,z)

所有的空间性质的计算器都包含在pywfn.gridprop子包下,其中每一个模块封装了一种类型的空间性质计算器

其包含的模块有

  • wfnfunc 分子轨道波函数
  • density 与电子密度相关的函数
  • dftgrid dft计算的空间格点

每个模块下都有一个Calculator类,实例化时传入要计算的分子即可

波函数

分子轨道波函数

计算指定空间格点指定分子轨道的分子轨道波函数,返回二维矩阵,形状为[nobt,npos](分子轨道数量,格点数量)

示例代码

下面代码是计算苯分子第21个分子轨道的20个随机点的波函数

python
from pywfn.base import Mole
from pywfn.gridprop import wfnfunc
import numpy as np

mole=Mole.from_file('./mols/C6H6.out')
caler=wfnfunc.Calculator(mole)

grids=np.random.rand(20,3)
val0,val1,val2=caler.obt_wfn(20,grids,2)
val0.shape,val1.shape,val2.shape
((20,), (20, 3), (20, 3, 3))

可以看到,我们计算了20个点,并且指定了计算level为2(这意味着需要计算一阶导和二阶导),计算得到三个数组,第一个数组是20个点的波函数的值,第二个数组是波函数一阶导,第二个是波函数二阶导

原子轨道波函数

计算指定空间格点的所有原子轨道的波函数

示例代码

下列代码计算苯分子空间中随机的10个点处102个原子轨道的波函数、波函数一阶导数和二阶导数

python
from pywfn.base import Mole
from pywfn.gridprop import wfnfunc
import numpy as np

mole=Mole.from_file("./mols/C6H6.out")
caler=wfnfunc.Calculator(mole)
grids=np.random.rand(10,3)
ato_wfns=caler.ato_wfn(grids,level=2)
print(len(ato_wfns)) # 102个原子轨道
ato_wfns[0][0].shape,ato_wfns[0][1].shape,ato_wfns[0][2].shape
102





((10,), (10, 3), (10, 3, 3))

可以看到,一共102个原子轨道,每个原子计算了10个点的原子轨道波函数数值、一阶导数和二阶导数

电子密度

分子总电子密度

计算空间中某些点处分子的总电子密度及电子密度的一阶导数和二阶导数

示例代码

python
from pywfn.base import Mole
from pywfn.gridprop import density
import numpy as np

mole=Mole.from_file("./mols/C6H6.out")
caler=density.Calculator(mole)

grids=np.random.rand(20,3)
dens0,dens1,dens2=vals=caler.mol_rho_dm(grids,level=2)

dens0.shape,dens1.shape,dens2.shape,
((20,), (20, 3), (20, 3, 3))

这里我们计算了随机20个点的电子密度,返回结果形式上与分子轨道波函数相同

预分子电子密度

将分子中原子的自由电子密度加和即为预分子电子密度,无需量化计算,只需要知道几何结构即可

示例代码

python
from pywfn.base import Mole
from pywfn.gridprop import density
import numpy as np

mole=Mole.from_file("./mols/C6H6.out")
caler=density.Calculator(mole)
grids=np.random.rand(20,3)
caler.pro_rho(grids)
[0.04626088986462919,
 0.050814897486494766,
 0.0503452744729373,
 0.08077468883028027,
 0.08410363597736738,
 0.06448774650780968,
 0.05613357776224836,
 0.0510525980079536,
 0.04582923193836833,
 0.04453544332498258,
 0.07550460275011656,
 0.06576113811026772,
 0.048480823198649596,
 0.03697888714168549,
 0.06348953798924459,
 0.05945919912735196,
 0.054704269572107764,
 0.036999873905360456,
 0.060066144063685066,
 0.03615585426662106]

RDG函数

约化密度梯度:Reduced Density Gradient

用来识别弱相互作用

计算公式

s=12(3π2)1/3|Δρ|ρ4/3

示例代码

python
from pywfn.base import Mole
from pywfn.gridprop import density
import numpy as np

mole=Mole.from_file("./mols/C6H6.out")
caler=density.Calculator(mole)
grids=np.random.rand(10,3)
caler.RDG(grids)
[100.0,
 0.7980327365008094,
 100.0,
 0.8274192655094126,
 100.0,
 0.818530917175119,
 100.0,
 100.0,
 100.0,
 0.8650417894032892]

其它的实空间函数的用法都是一样的(有的只能计算数值、有的能计算到二阶导数),在这里就不一一罗列了,可以在理论/实空间函数查看

DFT格点

dft计算的时候需要计算空间电子密度相关函数(泛函)的积分,如果对空间均匀采点会导致计算量非常大,所以采用非均匀采点

围绕每个原子采点,分为径向和球面采点两部分,有各自的权重,然后再把分子中所有原子的采点再归一化权重

首先展示一下径向采点

python
from pywfn.gridprop import dftgrid
dists,weits=dftgrid.rad_grids(1,20)
for dist,weit in zip(dists,weits):
    print(dist,weit)
178.06427461085983 142426688.49737817
44.017472640126535 1087888.3066196071
19.195669358089223 61292.43081822697
10.510047722828086 7748.70783895954
6.492092241699643 1512.53440507188
4.311941110422728 385.54559539720776
3.000000000000001 117.2205001597178
2.151298692346266 40.21441860480548
1.5724165284311622 14.993420940356227
1.1615314473505023 5.908529814152646
0.8609323512342596 2.405988775292347
0.635963805975586 0.9919622213710961
0.4648354984631968 0.4056743181625743
0.33333333333333354 0.1607962965153881
0.23191411347961657 0.05998446924111686
0.1540335477023657 0.02020213046411374
0.0951470465569794 0.005749117147523596
0.05209508360168708 0.0012251444061793398
0.02271825118574386 0.00014956627960667206
0.005615949646190353 4.468179459768279e-06

这里,dist是径向采点到原子的距离,可以看到离原子越近的地方采点越密集,但是相应地权重就越小

然后我们展示一下球面采点

python
from pywfn.gridprop import dftgrid
gridss,weits=dftgrid.sph_grids(6)
for grid,weit in zip(gridss,weits):
    print(grid,weit)
[1. 0. 0.] 0.16666666666666666
[-1.  0.  0.] 0.16666666666666666
[0. 1. 0.] 0.16666666666666666
[ 0. -1.  0.] 0.16666666666666666
[0. 0. 1.] 0.16666666666666666
[ 0.  0. -1.] 0.16666666666666666

sph_grids函数返回 Lebedev网格 ,这里返回了六个点的坐标及权重,格点的数量不是任意的,可选的值为6, 14, 26, 38, 50, 74, 86, 110, 146, 170, 194, 230, 266, 302, 350, 434, 590

有了径向格点和球面格点之后,我们就可以为每个原子一层一层的格点,然后组合成分子的dft格点

还是以苯环为例

python
from pywfn.base import Mole
from pywfn.gridprop import density,dftgrid
mole=Mole.from_file("./mols/C6H6.out")
grid_caler=dftgrid.Calculator(mole)
grids,weits=grid_caler.mol_grids(100,74)
dens_caler=density.Calculator(mole)
dens=dens_caler.mol_rho_dm(grids,level=0)[0]
np.sum(dens*weits)
np.float64(42.02238371780151)

这段代码我们生成分子格点,然后计算这些格点的电子密度,最后加权求和,实际上这就是对分子电子密度积分,得到的结果是42,苯环正是有42个电子