24小时热门版块排行榜    

北京石油化工学院2026年研究生招生接收调剂公告
查看: 3135  |  回复: 6
【奖励】 本帖被评价3次,作者beefly增加金币 2.4
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

[资源] 用数值方法计算偶极矩

此文目的:
1. 介绍偶极矩计算公式
2. 对于没有密度矩阵的理论方法,通过数值差分计算偶极矩

通常情况下,偶极矩都是用解析方法得到的。以CO分子为例,键长取1.2 Angstrom,用Gaussian程序计算CCSD/3-21G级别偶极矩的输入如下:

# ccsd/3-21g density=cc

Dipole moment with CCSD density

0 1
C
O 1 1.2


计算得到的偶极矩为-0.0743 Debye,方向为-Z(即,从O指向C)。很多初学者忘了加density=cc,这种情况下得到的是Hartree-Fock级别的偶极矩。

如果用数值方法,偶极矩该怎么算呢?偶极矩定义为能量对电场强度的一阶导数:(μx, μy, μz) = (dE/dFx, dE/dFy, dE/dFz),因此可以把偶极矩计算分解为不同电场下同一理论级别的一系列单点能计算。

从很多数值分析书上可以查到3点数值差分公式,
dE/dFx = [-E(-Dx) + E(+Dx)]/(2Dx)

其中Dx是X坐标轴方向上某个自定义的电场大小;E(-Dx)、E(+Dx)分别为给定电场下计算的分子总能量。还可以采用更精确的5点数值差分公式,
dE/dFx = [E(-2Dx) -8E(-Dx)+8E(+Dx)-E(+2Dx)]/(12Dx)

类似地,可以得到dE/dFy和dE/dFz的数值差分公式。

从上面的公式可以看到,3点公式只需要算两个点,5点公式只需要4个点。这是因为在一阶导数的数值差分计算中,用不到电场为0的点。

有了公式,下面可以设计计算流程了。仍是CO分子,C-O键仍放在Z轴。因为在X和Y方向,分子是对称的,偶极矩为零,故只需要算μz = dE/dFz。Gaussian提供了几种定义电场的方法,其中一种的命令为Field=z-n或Field=z+n,z表示电场方向,-n/+n表示电场为-/+n * 0.0001 a.u.。以下是第一个点的输入,电场取-0.002 a.u.:

# ccsd/3-21g Field=z-20

CCSD energy calculation with an electric dipole field in the Z direction of -0.002 au

0 1
C
O 1 1.2


第二、三、四个点的电场分别为-0.001,+0.001,+0.002 a.u.,需要分别把Field=z-20替换为Field=z-10,Field=z+10,和Field=z+20。四个点的CCSD能量如下(单位Hartree):

-112.2970354
-112.2970453
-112.2971037
-112.2971523


代入5点公式(Dz=0.001),换算成Debey(换算因子2.542)后,为-0.0742 Debey,非常接近上面通过解析方法得到的-0.0743 Debey。如果用三点公式(取第2、3个点),仍然得到-0.0742 Debey。

应用

上面计算偶极矩的数值方法有什么用呢?我们知道,Gaussian程序目前还不支持CCSD(T),BCCD(T),QCISD(TQ),等理论级别的密度矩阵。如果想得到这些理论级别的偶极矩,要么换程序(例如免费的CFour程序),要么就用上面的数值方法。

还有些理论方法是没有密度矩阵的,例如CBS,BSSE,Gn理论等。用上面的数值方法,原则上也可以计算这些理论级别的偶极矩。

[ Last edited by beefly on 2011-6-26 at 11:56 ]

补充:
全电子相对论方法(如ZORA,DKH2,等)的哈密顿需要做图景变换,因此一般情况下不能用上面的数值方法算偶极矩,除非在电场计算程序中加入图景变换。

[ Last edited by beefly on 2012-8-2 at 21:56 ]
回复此楼

» 收录本帖的淘帖专辑推荐

@个人收集@ 我的收藏

» 猜你喜欢

» 本主题相关商家推荐: (我也要在这里推广)

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
引用回帖:
6楼: Originally posted by 严登鑫 at 2017-08-03 11:21:16
好的,谢谢beefly~~

我需要得到电场中的最低能量,#p B3LYP/6-311G(3df) Field=x+60 Opt=Z-Matrix NoSymm freq 可以得到zero-point energies 这样做正确吗?...

1,确定电场方向和分子方位是否一致

2,g09加电场的情况下,dft没有解析hessian。频率计算很耗时建议用g16做
7楼2017-08-03 16:53:10
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 7 个回答

mei1988925

金虫 (小有名气)


★★★★★ 五星级,优秀推荐

学习了!谢谢beefly!
2楼2011-06-26 13:08:43
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lihb734

铁杆木虫 (职业作家)


★★★★★ 五星级,优秀推荐

感谢您的无私贡献
3楼2011-06-26 13:09:09
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

匿名

用户注销 (正式写手)


本帖仅楼主可见
4楼2017-08-02 16:47:22
已阅   申请QC强帖   回复此楼   编辑   查看我的主页
☆ 无星级 ★ 一星级 ★★★ 三星级 ★★★★★ 五星级
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 一志愿a区211,085601-307分求调剂 +10 党嘉豪 2026-03-31 20/1000 2026-03-31 16:57 by JourneyLucky
[考研] 江苏苏北高校诚邀调剂同学 +3 zzll406 2026-03-31 3/150 2026-03-31 16:54 by 及时行乐fan
[基金申请] 面上5B能上会吗? +8 redcom 2026-03-29 8/400 2026-03-31 15:53 by niuailing
[考研] 一志愿华东师范大学有机化学专业,初试351分,复试被刷求调剂! +8 真名有冰 2026-03-29 9/450 2026-03-31 15:40 by whf1979
[考研] 085601 329分调剂 +6 yzsa12 2026-03-31 6/300 2026-03-31 15:23 by yanflower7133
[考研] 085600,321分求调剂 +4 大馋小子 2026-03-31 4/200 2026-03-31 14:53 by 六翼TMR
[考研] 求调剂 +8 11ggg 2026-03-30 8/400 2026-03-31 13:56 by nanaliuyun
[考研] 354求调剂 +3 lxb598 2026-03-31 4/200 2026-03-31 13:42 by sophie2180
[考研] 求调剂 +4 研s. 2026-03-25 4/200 2026-03-31 10:05 by fmesaito
[考研] 0817化工学硕调剂 +7 努力上岸中! 2026-03-31 7/350 2026-03-31 09:58 by nalakaiqi
[考研] 一志愿中海洋320化学工程与技术学硕求调剂 +8 披星河 2026-03-30 8/400 2026-03-31 08:53 by lbsjt
[考研] 材料科学与工程求调剂 +10 深V宿舍吧 2026-03-29 10/500 2026-03-31 04:08 by fmesaito
[考研] 11408总分309,一志愿东南大学求调剂,不挑专业 +5 天赋带到THU 2026-03-29 6/300 2026-03-30 20:49 by dick_runner
[考研] 317分 一志愿南理工材料工程 本科湖工大 求调剂 +12 芋泥小铃铛 2026-03-28 12/600 2026-03-30 17:06 by wangjy2002
[考研] 275求调剂 +15 Micky11223 2026-03-25 20/1000 2026-03-29 20:44 by 唐沐儿
[考研] 299求调剂 +10 15188958825 2026-03-25 10/500 2026-03-29 17:51 by 王亮_大连医科大
[考研] 283求调剂 +7 A child 2026-03-28 7/350 2026-03-28 12:05 by zllcz
[考研] 0703化学求调剂,各位老师看看我!!! +5 祁祺祺 2026-03-25 5/250 2026-03-27 21:44 by 东方猪猪
[考研] 085602化学工程求调剂。 +4 平乐乐乐 2026-03-26 4/200 2026-03-26 17:57 by fmesaito
[考研] 303求调剂 +6 蓝山月 2026-03-25 6/300 2026-03-25 22:47 by 418490947
信息提示
请填处理意见