新手使用opt freq优化并计算频率，输出文件很久不更新也不报错

dravendr

近期老板让做一个体系，我需要用高斯算一下一个分子结晶以后的体系，算一下这种情况下的能量。我的方法是这样的:在剑桥数据库上下载下来了cif文件，然后用gview打开，保存成gjf后手动输入关键词。按照sob老师算例里最常见的计算级别B3LYP/6-31G*，先做一个结构优化和频率以求得热力学矫正量，然后换高级别算单点。我的输入文件如下：
%chk=DNTF crystal 283.15.chk
%mem=64GB
# B3LYP/6-31g* opt freq temperature=283.15

Title Card Required

0 1
(约80个C N O 原子)
Tv                 6.66200000    0.00000000    0.00000000
Tv                 0.00000000   10.74000000    0.00000000
Tv                 0.00000000    0.00000000   15.09300000

由于我除了上课，回来以后没有任何计算经验，所以不知道是否合理，请各位老师批评指正现在碰到的一个问题是输出文件6个多小时里，一直不更新，但是也没有报错，chk文件从最后写入的时间来看还是有在更新的，但是可以从截图（图1）里看到这个文件莫名其妙的非常大，我用formchk试着转换一下chk下下来放到gview里面看，但是他会给一个错误（图2）

out文件的最后是这么写的：

  1320 basis functions,  2464 primitive gaussians,  1320 cartesian basis functions
   312 alpha electrons      312 beta electrons
       nuclear repulsion energy     17236.6277922011 Hartrees.
NAtoms=   88 NActive=   88 NUniq=   88 SFac= 1.00D+00 NAtFMM=   60 NAOKFM=T Big=T
Integral buffers will be    131072 words long.
Raffenetti 2 integral format.
Two-electron integral symmetry is turned off.
FOutLm=  100.00.
Periodicity:                           1               1               1
Max integer dimensions:                8               5               4
   PBC vector 1 X=    12.5894 Y=    0.0000 Z=    0.0000
   PBC vector 2 X=     0.0000 Y=   20.2957 Z=    0.0000
   PBC vector 3 X=     0.0000 Y=    0.0000 Z=   28.5216
  Recp vector 1 X=     0.0794 Y=    0.0000 Z=    0.0000
  Recp vector 2 X=     0.0000 Y=    0.0493 Z=    0.0000
  Recp vector 3 X=     0.0000 Y=    0.0000 Z=    0.0351
Generated k point mesh (from -Pi to Pi):
K space mesh: X=      22 Y=      14 Z=      10
A half-cell shift: 0
Using k point mesh (from -Pi to Pi):
K space mesh: X=      22 Y=      14 Z=      10
A half-cell shift: 0
CountK=T Total number of k points:       0
CountK=T Total number of k points:    1544
One-electron integrals computed using PRISM.
NBasis=  1320 RedAO= T EigKep=  2.88D-04  NBF=  1320
NBsUse=  1320 1.00D-06 EigRej= -1.00D+00 NBFU=  1320
RepCel:  MaxNCR=     157 NClRep=     157 NMtPBC=    1267.
Not enough memory for dipole integrals, short by  1345199050 words.
Harris functional with IExCor=  402 and IRadAn=       5 diagonalized for initial guess.
RepCel:  MaxNCR=     157 NClRep=     157 NMtPBC=    1267.
HarFok:  IExCor=  402 AccDes= 0.00D+00 IRadAn=         5 IDoV= 2 UseB2=F ITyADJ=14
ICtDFT=  3500011 ScaDFX=  1.000000  1.000000  1.000000  1.000000
FoFCou: FMM=T IPFlag=      524288 FMFlag=      990000 FMFlg1=        3001
         NFxFlg=           0 DoJE=T BraDBF=F KetDBF=T FulRan=T
         wScrn=  0.000000 ICntrl=       500 IOpCl=  0 I1Cent=   200000004 NGrid=           0
         NMat0=    1 NMatS0=      1 NMatT0=    0 NMatD0=    1 NMtDS0=    0 NMtDT0=    0
Symmetry not used in FoFCou.
FMM levels:   4  Number of levels for PrismC:   3

-------------------- ATTENTION--------------------
Total charge is not ZERO ---      -0.000000000222
-------------------- ATTENTION--------------------

Requested convergence on RMS density matrix=1.00D-07 within 128 cycles.
Requested convergence on MAX density matrix=1.00D-05.
Requested convergence on             energy=1.00D-05.
No special actions if energy rises.
Diagonalized old Fock matrix for initial guess.
RepCel:  MaxNCR=     157 NClRep=     157 NMtPBC=    1267.

我在论坛里看到有前辈碰到类似的问题，说是计算频率会非常耗时，请问我目前的情况也是类似的吗，还是需要更改一下输入结构，恳请各位老师指点

sobereva

你的输入文件有Tv（平移矢量），Gaussian当成了周期性体系来计算，这么大体系根本算不动。Gaussian的周期性计算对于三维体系来说完全是个垃圾，没有任何实用性，速度超慢

你若要用Gaussian算，要么算单个分子，要么算多个分子的簇。
如果你要当周期性来算整个晶胞，推荐用CP2K
使用Multiwfn非常便利地创建CP2K程序的输入文件
http://sobereva.com/587（http://bbs.keinsci.com/thread-21668-1-1.html）

dravendr

sobereva 发表于 2021-11-8 21:55
你的输入文件有Tv（平移矢量），Gaussian当成了周期性体系来计算，这么大体系根本算不动。Gaussian的周期性 ...

感谢sob老师的解答，不过我的任务是想要对比这个分子在某个温度下，是容易结晶还是容易在一个混合溶剂里面溶解，所以我的思路是使用高斯来计算前者和后者。后者使用隐式溶剂模型（由于混合溶剂中使用的分子比较复杂，我查了一下高斯手册里面木有）或者是有机物混在一起的多个分子的簇。如果按照您说的方法，我使用cp2k或者vasp这种算周期性体系的软件的话，能量应该没法很好地相互对比吧。
所以我现在的思路是这样的，A,B,C三种分子都较大。对比分子簇的能量来评估A是否容易在B/C混合溶剂中结晶（B:C=1：1）：
对于结晶状态，我计算5个A分子的簇优化后的能量，得到的总值除以5，此体系记为A_cry；对于溶解状态，我计算1个A分子，两个B分子，两个C分子的簇优化过后的能量，此体系记为A2B2C_sol;然后我再计算两个B分子和两个C分子的簇优化过后的能量，此体系记为2B2C_sol；这三个体系，先opt freq，再把这个opt freq过后的结构改级别，算出高级别下的单点能
然后就是三个体系不同温度下的热力学矫正量使用shermo读取opt freq结果计算，加上高级别单点能，获得不同温度下三种簇的能量
然后就是对比EA_cry+E2B2C_sol与EA2B2C_sol在不同温度下的变化趋势，寻找交点来获得见解请老师给出一些建议，多谢

sobereva

dravendr 发表于 2021-11-9 09:39
感谢sob老师的解答，不过我的任务是想要对比这个分子在某个温度下，是容易结晶还是容易在一个混合溶剂里 ...

你这种问题最适合跑GROMACS分子动力学来研究

dravendr

sobereva 发表于 2021-11-10 01:57
你这种问题最适合跑GROMACS分子动力学来研究

非常感谢！下次老师开GROMACS班是什么时候？一定报名来学！

sobereva

dravendr 发表于 2021-11-10 10:57
非常感谢！下次老师开GROMACS班是什么时候？一定报名来学！

明年上旬，可以关注“北京科音”微信公众号，每次报名之前至少一个礼拜都会通知报名时间