gaussian和MOKIT联用脚本几何优化失败（已解决，脚本可自取） - 计算化学公社

#!/path/to/python python3
import sys
import os
import json
import re
from xml.dom.minidom import parse
# parse settings file
with open("gp.json") as fo:
settings = json.load(fo)
G_NPROC = settings["Gaussian"]["nproc"]
G_MEMORY = settings["Gaussian"]["memory"]
G_LEVEL = settings["Gaussian"]["level"]
G_BASIS = settings["Gaussian"]["basis"]
M_NPROC = settings["mokit"]["nproc"]
M_MEMORY = settings["mokit"]["memory"]
M_LEVEL = settings["mokit"]["level"]
M_COMMANDS = settings["mokit"]["commands"]
M_TOFIND = settings["mokit"]["tofind"]
if not os.path.exists("tmp"):
os.mkdir("tmp")
os.chdir("tmp")
# Gaussian UHF calculation input file templete
MOKIT_INP = f"""%nproc={G_NPROC}
%mem={G_MEMORY}
%chk=mokit.chk
{G_LEVEL}
MOKIT TASK
"""
ANG2BOHR = 1.8897259886
ELEMENTS = [
"Bq",
"H ", "He",
"Li", "Be", "B ", "C ", "N ", "O ", "F ", "Ne",
"Na", "Mg", "Al", "Si", "P ", "S ", "Cl", "Ar",
"K ", "Ca", "Sc", "Ti", "V ", "Cr", "Mn", "Fe", "Co", "Ni", "Cu", "Zn", "Ga", "Ge", "As", "Se", "Br", "Kr",
"Rb", "Sr", "Y ", "Zr", "Nb", "Mo", "Te", "Ru", "Rh", "Pd", "Ag", "Cd", "In", "Sn", "Sb", "Te", "I ", "Xe",
"Cs", "Ba", "La", "Ce", "Pr", "Nd", "Pm", "Sm", "Eu", "Gd", "Tb", "Dy", "Ho", "Er", "Tm", "Yb", "Lu", "Hf", "Ta", "W ", "Re", "Os", "Ir", "Pt", "Au", "Hg", "Tl", "Pb", "Bi", "Po", "At", "Rn",
"Fr", "Ra", "Ac", "Th", "Pa", "U ", "Np", "Pu", "Am", "Cm", "Bk", "Cf", "Es", "Fm", "Md", "No", "Lr", "Rf", "Db", "Sg", "Bh", "Hs", "Mt", "Ds", "Rg", "Cn", "Nh", "Fl", "Mc", "Lv", "Ts", "Og",
]
(layer, InputFile, OutputFile, MsgFile, FChkFile, MatElFile) = sys.argv[1:]
atoms = []
coords = []
with open(InputFile, "r") as fi:
(atoms, derivs, charge, spin) = [int(x) for x in fi.readline().split()]
MOKIT_INP += f"{charge} {spin}\n"
for i in range(0, atoms):
arr = fi.readline().split()
atom = ELEMENTS[int(arr[0])]
coord = [f"{(float(x) / ANG2BOHR)}" for x in arr[1:4]]
MOKIT_INP += f" {atom}"
MOKIT_INP += " " * 4
MOKIT_INP += (" " * 4).join(coord)
MOKIT_INP += "\n"
MOKIT_INP += "\n"
MOKIT_INP = f"""{G_BASIS}
"""
with open("MOKIT.gjf", "w") as fo:
fo.write(MOKIT_INP)
print(">>> Starting UHF Calculation...")
os.system("g16 MOKIT.gjf")
os.system("formchk MOKIT.chk")
os.system("mv MOKIT.fchk MOKIT.fch")
print(">>> UHF Calculation Done!")
print(">>> Preparing automr input file...")
NEVPT2_INP = f"""%nproc={M_NPROC}
%mem={M_MEMORY}
{M_LEVEL}
mokit{{ {M_COMMANDS} }}
0 1
C 0.0 0.0 0.0
"""
with open("NEVPT2.gjf", "w") as fo:
fo.write(NEVPT2_INP)
# call MOKIT to do automr calculation
print(">>> Starting automr Calculation...")
os.system(f"automr NEVPT2.gjf >& NEVPT2.out")
print(">>> automr Calculation Done!")

复制代码

{
"Gaussian": {
"nproc": "8",
"memory": "32GB",
"level": "#p uhf gen stable=opt scf=(noincfock,novaracc,xqc) nosymm int=(nobasistransform)",
"basis":[
"n o",
"6-31g*",
"****",
"",
""
]
},
"mokit": {
"nproc": "8",
"memory": "32GB",
"level":"casscf(5,5)/gen",
"commands": "casscf_prog=pyscf,force,ist=2,readuhf="mokit.fch"",
"tofind": "E(CASSCF) ="
}
}

复制代码

#p CASPT2/def2TZVP
mokit{force}
电荷自旋多重度
坐标

复制代码

#p CASPT2/def2TZVP
mokit{force}
电荷自旋多重度
坐标

复制代码

Computing numerical gradient using default procedure for command NEVPT2
*** Long output written to logfile /szp/soft/scripts/gm/test2/tmp/UHF_uno2CASSCF_NEVPT2.log ***
Running default procedure: CASSCF000 NEVPT2
? Error
? Numerical calculation of lattice gradient impossible
? The problem occurs in driver(forces)
GLOBAL ERROR fehler on processor 0

复制代码