高斯计算(Gaussian)Input File模板

使用Gaussian软件的量子化学计算流程首先先简略看一遍用Gaussian做计算的流程。本文内容称不上是业内宝典名师权威，但对于那些没有充分了解量子化学但还想做量子化学计算的同学们来说是个很好的指导教程！我观察了我周围的情况后发现，很多委托别人来做Gaussian计算的人都是出于担心使用软件需要数学知识和计算能力，自己不能胜任的心理所以选择交给其他人来做。当然理解计算和数学原理是最好的，但不妨先姑且学学怎么用这个软件。Gaussian的计算流程十分简单：1.首先创建分子的结构信息（.xyz﹑笛卡尔坐标或z-matrix），可以在GaussView,Mercury,Avogadoro,Chem3D等软件上画出分子的结构。2.创建GaussianInputfile(.gjf)（在textedit或GaussView上操作）3.向Gaussian提交GaussianInputfile（在Gaussian上操作）4.验收・分析GaussianOutputfile（文件后缀名为.chk/.fchk/.log），可以在textedit,GaussView,iqmol上操作本文的全部内容都是关于整套流程中最重要的第二步——创建GaussianInputfile(.gjf)。GaussianInputfile实质上是一种文本文档，它汇总了我们想让Gaussian做什么样的计算的指令信息，这些指令信息在Inputfile上为软件可以读取的形式。我曾幻想过有一天某位神人能打造出一个完美的Inputfile，完美到让Gaussian直接解出薛定谔方程，那量子化学界也就迎来了它生命的终点，当然这是不可能的。那么，Inputfile该怎么写呢？通过在GaussView的操作，Inputfile可以自动生成，想要完美的Inputfile可以打开Inputfile的Text形式，在Text上修改。虽说规则都在HPC手册和Gaussian官方网站。但没有示例很难自己写可怎么办呢？不如到这篇文章里找找模板吧。只需将本文里的模板直接复制粘贴到textfile里面就好（因为格式差异，可能有些空行的长短有错误，还请读者见谅）。修改分子的XYZ轴的信息和最开始的两行，最后加上两行的空行差不多就可以投入到Gaussian里进行计算了。例1)以b3lyp/6-31g(d,p)为泛函/基组的基态计算（S0构造最优化，振动数计算）%nprocshared=28

%mem=100GB

%chk=TNT.chk

#optfreqb3lyp/6-31g(d,p)****01

C-4.592550.08304-0.05852

C-3.073180.17974-0.06585

C-2.358221.43873-0.06606

N-2.970162.74777-0.05954

O-2.216373.82816-0.06039

O-4.149592.90090-0.05332

C-0.948371.40993-0.07316

C-0.220550.22236-0.07994

N1.209960.24774-0.08700

O1.859231.39310-0.08718

O1.83665-0.76914-0.09283

C-0.91178-0.98707-0.07972

C-2.31693-1.04992-0.07289

N-2.91158-2.36710-0.07344

O-2.14448-3.43785-0.08011

O-4.09243-2.53018-0.06800

H-5.188630.99767-0.05309

H-4.91946-0.459090.85345

H-4.92857-0.45422-0.97006

H-0.394162.33820-0.07345

H-0.33521-1.90221-0.08503下面我们来看看上面的这个例子里每句话说的是什么意思。%nprocshared=28（使用的CPU个数，这个是根据你用的电脑决定的）

%mem=100GB（使用的内存值，这也是根据你用的电脑决定的）

%chk=TNT.chk（作成的checkpointfile即输出文件的名字）

#optfreqb3lyp/6-31g(d,p)（routesection从这一行起开始，这一串字符声明了你的计算条件。opt=结构最优化，freq＝振动数以b3lyp/6-31g(d,p)为泛函/基组进行计算）（这里的空行代表

routesection的结束）****（这一行是commentsection你可以输入你将要计算的分子的名字之类的）（这里的空行代表commentsection的结束）01（代表电荷和自旋多重度。两个数字之间一定要空一格。01代表电荷为0，多重度是1）

C-4.592550.08304-0.05852（下面的一大串数字表示分子中各个原子的XYZ轴的坐标信息，也有的inputfile是用Z-matrix来表示原子位置的）

C-3.073180.17974-0.06585

C-2.358221.43873-0.06606

N-2.970162.74777-0.05954

O-2.216373.82816-0.06039

O-4.149592.90090-0.05332

C-0.948371.40993-0.07316

C-0.220550.22236-0.07994

N1.209960.24774-0.08700

O1.859231.39310-0.08718

O1.83665-0.76914-0.09283

C-0.91178-0.98707-0.07972

C-2.31693-1.04992-0.07289

N-2.91158-2.36710-0.07344

O-2.14448-3.43785-0.08011

O-4.09243-2.53018-0.06800

H-5.188630.99767-0.05309

H-4.91946-0.459090.85345

H-4.92857-0.45422-0.97006

H-0.394162.33820-0.07345

H-0.33521-1.90221-0.08503（这里的空行代表上面的分子结构信息输入完毕）（如果有需要增加的指令，在这一行输入，具体参照例8）（这里的空行表示上行的增加指令已经结束）（这里的空行表示inputfile到此全部结束）这就是一个inputfile里那些数字字母所要传达的信息。例2)以M062X/6-31g(d,p)为泛函/基组的激发态计算（S1构造最优化，振动数计算）%mem=100GB

%nprocshared=28

%chk=TNT.chk

#optfreq

td=(singlets,nstates=10,root=1)

b3lyp/6-31g(d,p)****01

C-4.592550.08304-0.05852

C-3.073180.17974-0.06585

C-2.358221.43873-0.06606

N-2.970162.74777-0.05954

O-2.216373.82816-0.06039

O-4.149592.90090-0.05332

C-0.948371.40993-0.07316

C-0.220550.22236-0.07994

N1.209960.24774-0.08700

O1.859231.39310-0.08718

O1.83665-0.76914-0.09283

C-0.91178-0.98707-0.07972

C-2.31693-1.04992-0.07289

N-2.91158-2.36710-0.07344

O-2.14448-3.43785-0.08011

O-4.09243-2.53018-0.06800

H-5.188630.99767-0.05309

H-4.91946-0.459090.85345

H-4.92857-0.45422-0.97006

H-0.394162.33820-0.07345

H-0.33521-1.90221-0.08503例3)以cam-B3LYP/6-31g(d,p)为泛函/基组的激发态计算（S0构造）%mem=100GB

%nprocshared=28

%chk=TNT.chk

#

sptd=(singlets,nstates=10,root=1)

cam-b3lyp/6-31g(d,p)

Guess=ReadGeom=AllCheckpoint****这里的“Guess=ReadGeom=AllCheckpoint”是个很好用的操作，可以将之前执行的计算结果和构造不做改动直接用到其他计算任务里。具体应用场景见下面的例4。例4)NTO计算（S1跃迁轨道的可视化）%mem=100GB

%nprocshared=28

%chk=TNT.chk

#M062X/6-31g(d,p)

Geom=AllCheckGuess=(Read,Only)Density=(Check,Transition=1)Pop=(Minimal,NTO,SaveNTO)01这是标准的模板。如果不想覆盖原有的chkfile的话，改成下面这样：%mem=100GB

%nprocshared=28%oldchk=TNT.chk

%chk=TNT_NTO.chk#M062X/6-31g(d,p)Geom=AllCheckGuess=(Read,Only)Density=(Check,Transition=1)Pop=(Minimal,NTO,SaveNTO)01这样，计算结果就不会覆盖原有的chkfile而是新输出一个chkfile。例5)以M062X/6-31g(d,p)为泛函/基组的激发态计算（T1结构最优化，振动数计算）%mem=100GB

%nprocshared=28

%chk=TNT.chk

#optfreq

td=(50-50,nstates=10,root=1)

b3lyp/6-31g(d,p)****01

C-4.592550.08304-0.05852

C-3.073180.17974-0.06585

C-2.358221.43873-0.06606

N-2.970162.74777-0.05954

O-2.216373.82816-0.06039

O-4.149592.90090-0.05332

C-0.948371.40993-0.07316

C-0.220550.22236-0.07994

N1.209960.24774-0.08700

O1.859231.39310-0.08718

O1.83665-0.76914-0.09283

C-0.91178-0.98707-0.07972

C-2.31693-1.04992-0.07289

N-2.91158-2.36710-0.07344

O-2.14448-3.43785-0.08011

O-4.09243-2.53018-0.06800

H-5.188630.99767-0.05309

H-4.91946-0.459090.85345

H-4.92857-0.45422-0.97006

H-0.394162.33820-0.07345

H-0.33521-1.90221-0.08503这不是针对T1态的构造优化，而是分子会趋向于落入能量最低的激发状态。通常情况下，分子会趋向于处于T1态。例5-2)不用td的T1态构造最优化%mem=100GB

%nprocshared=28

%chk=TNT.chk

#optfreqb3lyp/6-31g(d,p)****03

C-4.592550.08304-0.05852

C-3.073180.17974-0.06585

C-2.358221.43873-0.06606

N-2.970162.74777-0.05954

O-2.216373.82816-0.06039

O-4.149592.90090-0.05332

C-0.948371.40993-0.07316

C-0.220550.22236-0.07994

N1.209960.24774-0.08700

O1.859231.39310-0.08718

O1.83665-0.76914-0.09283

C-0.91178-0.98707-0.07972

C-2.31693-1.04992-0.07289

N-2.91158-2.36710-0.07344

O-2.14448-3.43785-0.08011

O-4.09243-2.53018-0.06800

H-5.188630.99767-0.05309

H-4.91946-0.459090.85345

H-4.92857-0.45422-0.97006

H-0.394162.33820-0.07345

H-0.33521-1.90221-0.08503注意，与S0的结构最优化和振动数计算的模板不同的地方在于粗体字换成了03。例6)只优化氢原子的构造（用于XRD结晶构造解析后的构造最优化）%nprocshared=28

%mem=100GB

%chk=TNT.chk

#optfreqb3lyp/6-31g(d,p)****01

C

-1

-4.592550.08304-0.05852

C

-1

-3.073180.17974-0.06585

C

-1

-2.358221.43873-0.06606

N

-1

-2.970162.74777-0.05954

O

-1

-2.216373.82816-0.06039

O

-1

-4.149592.90090-0.05332

C

-1

-0.948371.40993-0.07316

C

-1

-0.220550.22236-0.07994

N

-1

1.209960.24774-0.08700

O

-1

1.859231.39310-0.08718

O

-1

1.83665-0.76914-0.09283

C

-1

-0.91178-0.98707-0.07972

C

-1

-2.31693-1.04992-0.07289

N

-1

-2.91158-2.36710-0.07344

O

-1

-2.14448-3.43785-0.08011

O

-1

-4.09243-2.53018-0.06800

H

0

-5.188630.99767-0.05309

H

0

-4.91946-0.459090.85345

H

0

-4.92857-0.45422-0.97006

H

0

-0.394162.33820-0.07345

H

0

-0.33521-1.90221-0.08503原子右侧的数字-1声明固定的意思，0声明可以移动的意思。即上面的一大串程序表明计算机在对分子作结构最优化时只被允许移动H原子。如果需要计算的分子太庞大以至于画出来很费时间，我会用下面的Python代码准备gjffile。filename_input=“./A.gjf”

filename_output=“./output.gjf”start_line_num=9

output_text=“”withopen(filename_input,“r”)asf:

lines=f.readlines()forline_num,lineinenumerate(lines):

line=line.strip()#removebreaklineifline_num+1<start_line_num:#line_numstartsfrom0

output_text+=“{}\n”.format(line)#addbreaklineattheend

else:

line_array=[eforeinline.split(”“)ife!=“”]#splitinputlinebyspaceandmakeanarray

iflen(line_array)>1:

ifline_array[0]==“H”:

line_array.insert(1,str(0))#add“0”incaseof“H”

else:

line_array.insert(1,str(-1))#otherwise,add“1”line_output=“\t”.join(line_array)

output_text+=“{}\n”.format(line_output)#addbreaklineattheend#print(“===outputtext====”)

#print(output_text)withopen(filename_output,“w”)asf_out:

f_out.write(output_text)f=open(‘output.gjf’,‘a’,encoding=’utf-8′,newline=’\n’)

f.write(‘\n’)

f.write(‘\n’)

f.close()例7)丁烷的旋转﹑势能面的计算%chk=butane_scan2

%mem=100GB

%nproc=36

#popt=modredundantb3lyp/6-31g(d)pop=fullmemo01

C-1.81226530-0.33960474-0.18951040

H-2.75438476-0.79975943-0.40297839

H-1.10701678-1.089529830.10227189

C-1.298923080.38635154-1.44691537

H-2.003260750.26246181-2.24281724

H-0.35589216-0.02730756-1.73756778

C-1.129120841.88597588-1.14049402

H-0.425063522.00984515-0.34434097

H-0.772093122.39031163-2.01403661

H-2.072226532.29971720-0.85020142

C-1.984691040.676048870.95518044

H-2.689607231.426220550.66322903

H-2.341785640.171724691.82870237

H-1.042482461.135882491.16894671D11147S3610.000000例8)根据原子指定基函数时%nprocshared=8

%mem=8GB

%c