24小时热门版块排行榜    

查看: 3534  |  回复: 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的回帖

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强帖   回复此楼   编辑   查看我的主页
引用回帖:
4楼: Originally posted by 严登鑫 at 2017-08-02 16:47:22
您好请教一个问题,我想用gaussian计算外加电场对分子能量的影响。
具体步骤如下:
                              1.提交Gaussian 程序进行分子优化及频率计算,得到该分子的最稳定结构
         
             ...

1,步骤1没问题,步骤2未必需要opt,看你要的是电场中的单点能,还是电场中的最低能量

2,有用的能量值有好几个,真空中绝对零度下的电子总能量,前者的零点能校正,以及考虑环境影响的吉布斯自由能

3,频率计算才有zpe

» 本帖已获得的红花(最新10朵)

5楼2017-08-03 09:13:56
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

匿名

用户注销 (正式写手)


送红花一朵
本帖仅楼主可见
6楼2017-08-03 11:21:16
已阅   申请QC强帖   回复此楼   编辑   查看我的主页
引用回帖:
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的回帖
相关版块跳转 我要订阅楼主 beefly 的主题更新
☆ 无星级 ★ 一星级 ★★★ 三星级 ★★★★★ 五星级
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[找工作] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 EZ67DG6CvcKL 2026-08-31 4/200 2026-09-01 06:06 by o3EMxHUAQDtS
[基金申请] 怎么看青基中了没有啊 +6 叶九微 2026-08-26 6/300 2026-08-31 23:54 by yudaoqian88
[基金申请] 国社科又开始会评了,不知道这次命运如何 +7 雨打竹帘 2026-08-30 11/550 2026-08-31 23:16 by hittle2008
[基金申请] 投票:  有多少人是今天查系统知道结果的? +17 爱看书的可乐 2026-08-26 19/950 2026-08-31 21:30 by xiangy672
[基金申请] 基金未中,这种答复是模板吗? +6 zhaosm1982 2026-08-27 7/350 2026-08-31 21:18 by qdxxmc
[基金申请] 面上意见出来了 +12 黄鸟于飞Chao 2026-08-29 23/1150 2026-08-31 18:57 by 黄鸟于飞Chao
[基金申请] 中青基了要发朋友圈吗? +7 349506619 2026-08-28 7/350 2026-08-31 13:39 by 冼亮淀粉酶
[基金申请] 29号明天会评吗 +4 笨笨唐 2026-08-28 4/200 2026-08-31 09:30 by huixian257
[基金申请] 我就是申请一个面上项目而已,这评审意见是按照杰青的条件评的吧? +6 gouxfjh 2026-08-28 11/550 2026-08-30 07:57 by gouxfjh
[基金申请] 国自然面上复盘~欢迎讨论 (金币+15) +15 晴天加油 2026-08-26 16/800 2026-08-29 18:28 by symmetry
[基金申请] 国自然评审意见 +13 wangmingqi 2026-08-28 19/950 2026-08-29 10:22 by Poppy1104
[基金申请] 2026年叶企孙基金 +4 bud_bud 2026-08-27 7/350 2026-08-29 07:23 by foolishmani
[基金申请] 系统查不到 +11 董八千 2026-08-26 11/550 2026-08-28 18:06 by Leogzhya
[基金申请] 基金系统什么内容也没有 30+4 winsaint 2026-08-27 9/450 2026-08-28 11:06 by maolC
[基金申请] 看板上这么多中的,有点像50人群里49个人都是骗子的那种感觉…… +5 a089 2026-08-26 6/300 2026-08-27 14:05 by jonewore
[基金申请] 我不理解! +15 Edward_pc 2026-08-26 23/1150 2026-08-26 20:34 by zzuzxg
[基金申请] 为什么国自然不能直接公布 +4 bjdxyxy 2026-08-26 4/200 2026-08-26 13:12 by qingmu1201
[基金申请] 系统进不去 +4 yanglien 2026-08-26 5/250 2026-08-26 11:10 by wenfengw83
[基金申请] 今天务委会开完了,明天出结果吗 +19 angus9576 2026-08-25 23/1150 2026-08-26 10:03 by zp519
[基金申请] 在坚冰还盖着北海的时候,我看到了怒放的梅花。 (金币+10) +6 ziyangfang 2026-08-25 9/450 2026-08-25 20:26 by huagongfeihu
信息提示
请填处理意见