6.20. 复合方法 ¶
复合方法是一种复杂的脚本语言,可以直接在 ORCA 的输入中使用。使用“复合”,用户可以将正常 ORCA 计算的各个部分组合起来,以评估自定义函数。为了说明其用法,我们将使用一个示例。有关该模块的更详细描述,请参阅复合方法部分。
6.20.1. 示例
作为一个典型的例子,我们将使用描述
# ----------------------------------------------
# Umbrella coordinate mapping for NH3
# Author: Frank Neese
# ----------------------------------------------
variable JobName = "NH3-umbrella";
variable amin = 50.0;
variable amax = 130.0;
variable nsteps = 21;
Variable energies[21];
Variable angle;
Variable JobStep;
Variable JobStep_m;
variable step;
Variable method = "BP86";
Variable basis = "def2-SVP def2/J";
step = 1.0*(amax-amin)/(nsteps-1);
# Loop over the number of steps
# ----------------------------
for iang from 0 to nsteps-1 do
angle = amin + iang*step;
JobStep = iang+1;
JobStep_m= JobStep-1;
if (iang>0) then
Read_Geom(JobStep_m);
New_step
! &{method} &{basis} TightSCF Opt
%base "&{JobName}.step&{JobStep}"
%geom constraints
{A 1 0 2 &{angle} C}
{A 1 0 3 &{angle} C}
{A 1 0 4 &{angle} C}
end
end
Step_End
else
New_step
! &{method} &{basis} TightSCF Opt
%base "&{JobName}.step&{JobStep}"
%geom constraints
{A 1 0 2 &{angle} C}
{A 1 0 3 &{angle} C}
{A 1 0 4 &{angle} C}
end
end
* int 0 1
N 0 0 0 0.0 0.0 0.0
DA 1 0 0 2.0 0.0 0.0
H 1 2 0 1.06 &{angle} 0.0
H 1 2 3 1.06 &{angle} 120.0
H 1 2 3 1.06 &{angle} 240.0
*
Step_End
endif
Read energies[iang] = SCF_ENERGY[jobStep];
print(" index: %3d Angle %6.2lf Energy: %16.12lf Eh\n", iang, angle, energies[iang]);
EndFor
# Print a summary at the end of the calculation
# ---------------------------------------------
print("////////////////////////////////////////////////////////\n");
print("// POTENTIAL ENERGY RESULT\n");
print("////////////////////////////////////////////////////////\n");
variable minimum,maximum;
variable Em,E0,Ep;
variable i0,im,ip;
for iang from 0 to nsteps-1 do
angle = amin + 1.0*iang*step;
JobStep = iang+1;
minimum = 0;
maximum = 0;
i0 = iang;
im = iang-1;
ip = iang+1;
E0 = energies[i0];
Em = E0;
Ep = E0;
if (iang>0 and iang<nsteps-1) then
Em = energies[im];
Ep = energies[ip];
endif
if (E0<Em and E0<Ep) then minimum=1; endif
if (E0>Em and E0>Ep) then maximum=1; endif
if (minimum = 1 ) then
print(" %3d %6.2lf %16.12lf (-)\n",JobStep,angle, E0 );
endif
if (maximum = 1 ) then
print(" %3d %6.2lf %16.12lf (+)\n",JobStep,angle, E0 );
endif
if (minimum=0 and maximum=0) then
print(" %3d %6.2lf %16.12lf \n",JobStep,angle, E0 );
endif
endfor
print("////////////////////////////////////////////////////////\n");
End # end of compound block
让我们先看看如何执行这个输入。为了运行它,最简单的方法是将其保存在一个普通文本文件中,使用名称“umbrella.cmp”,然后使用以下 ORCA 输入文件:
%Compound "umbrella.cmp"
不需要更多。ORCA 将读取复合文件并采取适当的行动。
关于这个 ORCA 输入的几点说明。首先,没有简单的输入行(以“!”开头)。使用复合功能时不需要简单输入,但如果用户添加了简单输入,所有来自简单输入的信息将传递给实际的复合作业。
此外,如果不想创建单独的复合文本文件,在 ORCA 中完全可以像使用其他 ORCA 块一样使用复合功能。这意味着在%Compound 指令之后,可以附加复合文件的内容,而不是给出文件名。
正如我们将看到的,在复合脚本文件中,每个复合作业可以包含正常 ORCA 输入文件的所有信息。这里有两个非常重要的例外:处理器的数量和 MaxCore。这些信息应在初始 ORCA 输入文件中设置,而不是在实际的复合文件中。
复合块的结构与所有 ORCA 块相同。它以“%”开头,以“End”结尾,如果输入不是从文件中读取的。如果复合指令在文件中,如上面的例子所示,则只需在括号内提供文件名,而不需要最后的 END。
6.20.2. 定义变量 ¶
正如我们已经指出的,可以在复合块内提供所有计算和数据处理的信息,或者创建一个包含所有细节的普通文本文件,让 ORCA 读取。后者的优点是可以使用同一个文件处理多个几何体。在前面的例子中,我们将 ORCA 指向一个外部文件。该文件“umbrella.cmp”包含所有必要的信息。
现在让我们尝试分析复合“umbrella.cmp”文件。
# ----------------------------------------------
# Umbrella coordinate mapping for NH3
# Author: Frank Neese
# ----------------------------------------------
variable JobName = "NH3-umbrella";
variable amin = 50.0;
variable amax = 130.0;
variable nsteps = 21;
Variable energies[21];
Variable angle;
Variable JobStep;
Variable JobStep_m;
variable step;
Variable method = "BP86";
Variable basis = "def2-SVP def2/J";
step = 1.0*(amax-amin)/(nsteps-1);
第一部分包含一些一般性评论和变量定义。对于评论,我们使用与普通 ORCA 输入相同的语法,通过“#”符号。请注意,同一行中多个“#”符号会导致错误。
在初始评论之后,我们看到了一些声明和定义。声明变量的方式有很多种,详细描述见“变量 - 一般”部分。
所有变量声明以指令 Variable 开头,这表示程序将期待一个或多个新变量的声明。然后有许多选项,包括定义多个变量、为变量赋值或使用值列表。然而,所有声明必须以;符号结束。这个符号是程序的一个消息,表示当前命令的结束。在 Compound 中,每个命令末尾需要;符号是一个普遍要求,只有极少数例外。
6.20.3. 运行计算 ¶
# Loop over the number of steps
# ----------------------------
for iang from 0 to nsteps-1 do
angle = amin + iang*step;
JobStep = iang+1;
JobStep_m= JobStep-1;
if (iang>0) then
Read_Geom(JobStep_m);
New_step
! &{method} &{basis} TightSCF Opt
%base "&{JobName}.step&{JobStep}"
%geom
constraints
{A 1 0 2 &{angle} C}
{A 1 0 3 &{angle} C}
{A 1 0 4 &{angle} C}
end
end
Step_End
else
New_step
! &{method} &{basis} TightSCF Opt
%base "&{JobName}.step&{JobStep}"
%geom
constraints
{A 1 0 2 &{angle} C}
{A 1 0 3 &{angle} C}
{A 1 0 4 &{angle} C}
end
end
* int 0 1
N 0 0 0 0.0 0.0 0.0
DA 1 0 0 2.0 0.0 0.0
H 1 2 0 1.06 &{angle} 0.0
H 1 2 3 1.06 &{angle} 120.0
H 1 2 3 1.06 &{angle} 240.0
*
Step_End
endif
Read energies[iang] = SCF_ENERGY[jobStep];
print(" index: %3d Angle %6.2lf Energy: %16.12lf Eh\textbackslash{}n", iang, angle, energies[iang]);
EndFor
然后我们进入信息最密集的部分。我们从 for 循环的定义开始。复合 for 循环的语法是:
对于变量 From startValue 到 endValue 执行
指令
EndFor
如上例所示,startValue 和 endValue 可以是常数或之前定义的变量,甚至是这些变量的函数。请记住,它们必须是整数。循环结束的信号是 EndFor 指令。有关 for 循环的更多细节,请参阅 For 部分。
然后我们继续分配一些变量。
angle = amin + iang*step;
JobStep = iang+1;
JobStep_m = JobStep-1;
变量赋值的语法与每种编程语言中的变量相似,后面跟着 = 符号,然后是值或方程。请记住,赋值必须始终以 ; 符号结束。
下一步是每种编程语言中另一个重要部分,即 if 块。if 块的语法如下:
如果 (要评估的表达式) 那么
指令
否则如果 ( 要评估的表达式 ) 那么
指令
否则
指令
EndIf
else if 和 else 部分是可选的,但最终的 EndIf 必须始终标志着 if 块的结束。有关 if 块使用的更多细节,请参阅手册的 If 部分。
接下来我们有一个特定于复合方法的命令,而不是普通编程语言的一部分。这是 ReadGeom 命令。它的语法是:
Read_Geom(整数值);
在解释这个命令之前,我们将继续处理复合脚本中的下一个命令,然后再返回到这个命令。
下一个命令是所有复合脚本的基础。这是 New_Step 命令。该命令向复合指示接下来是一个正常的 ORCA 计算。它的语法是:
New_Command 正常 ORCA 输入 Step_End
关于 New_Step 命令的一些评论。首先,在 New_Step - Step_End 命令内部,可以添加所有正常 ORCA 输入接受的可能命令。我们应该记住,定义处理器数量的命令和 MaxCore 命令将被忽略。
另一个需要记住的要点是步骤的概念。每个 New_Step - Step_End 结构对应一个步骤,从 1 开始计数(第一个 ORCA 计算)。这有助于我们定义该计算将创建的属性文件,以便我们可以用它来检索信息。
在 New_Step - Step_End 块中的一个重要特性是使用结构 &{variable}。该结构允许用户在 New_Step - Step_End 块内部使用在外部定义的变量,从而使 ORCA 输入更加通用。例如,在上面的脚本中,我们构建了计算的基本名称。
%base "&{JobName}.step&{JobStep}"
使用定义的变量 JobName 和 JobStep。有关 &{} 结构的使用的更多细节,请参阅 & 部分,而有关 New_Step - Step_End 结构,请参阅 New_Step 部分。
最后,关于计算的几何形状有几点评论。提供几何形状给 New_Step - Step_End 计算有三种方式。第一种是传统的 ORCA 输入方式,我们可以像在所有 ORCA 输入中那样提供坐标或包含坐标的文件名。然而在 Compound 中,如果我们没有传递任何关于计算几何形状的信息,那么 Compound 将自动尝试读取上一个计算的几何形状。这是提供几何形状给复合步骤的第二种(隐式)方式。然后还有第三种方式,这就是我们在上面的例子中使用的。这是 Read_Geom 命令。该命令的语法是:
Read_Geom (步骤编号);
当我们想要将特定的几何体传递给一个计算时,可以使用此命令,该几何体并未在 New_Step - Step_End 结构中明确给出,也不是来自前一步。然后,我们只需在运行新计算之前传递我们感兴趣的计算步骤的编号。有关 Read_Geom 命令的更多详细信息,请参阅 Read_Geom 部分。
6.20.4. 数据处理 ¶
Compound 的一个最强大的功能是它可以直接访问计算的属性。为了使用这些属性,我们定义了 Read 命令。在前面的例子中,我们用它来读取计算的 SCF 能量:
Read energies[iang] = SCF\_ENERGY[jobStep];
该命令的语法是:
读取变量名 = property
其中变量名称是已经定义的变量的名称,属性是来自已知属性表中的已知属性,步骤是我们感兴趣的计算步骤。有关 Read 命令的更多详细信息,请参阅 Read 部分。
# Print a summary at the end of the calculation
# ---------------------------------------------
print("////////////////////////////////////////////////////////\\n");
print("// POTENTIAL ENERGY RESULT\\n");
print("////////////////////////////////////////////////////////\\n");
variable minimum,maximum;
variable Em,E0,Ep;
variable i0,im,ip;
for iang from 0 to nsteps-1 do
angle = amin + 1.0*iang*step;
JobStep = iang+1;
minimum = 0;
maximum = 0;
i0 = iang;
im = iang-1;
ip = iang+1;
E0 = energies[i0];
Em = E0;
Ep = E0;
if (iang>0 and iang<nsteps-1) then
Em = energies[im];
Ep = energies[ip];
endif
if (E0<Em and E0<Ep) then minimum=1; endif
if (E0>Em and E0>Ep) then maximum=1; endif
if (minimum = 1 ) then
print(" %3d %6.2lf %16.12lf (-)\textbackslash{}n",JobStep,angle, E0 );
endif
if (maximum = 1 ) then
print(" %3d %6.2lf %16.12lf (+)\textbackslash{}n",JobStep,angle, E0 );
endif
if (minimum=0 and maximum=0) then
print(" %3d %6.2lf %16.12lf \textbackslash{}n",JobStep,angle, E0 );
endif
endfor
print("////////////////////////////////////////////////////////\\n");
一旦所有数据可用,我们可以像在任何编程语言中一样在方程中使用它们。
print 语句的语法是:
print( format string, [variables]);
例如在之前的脚本中,我们这样使用它:
print(" %3d %6.2lf %16.12lf \n",JobStep,angle, E0 );
其中 %3d、%6.2lf 和 %16.2lf 是格式标识符,JobStep、angle 和 E0 是之前定义的变量。语法与编程语言 C 中广泛接受的 printf 命令的语法非常相似。有关 print 语句的更多详细信息,请参阅章节:Print。
与 print 命令类似的是 write2file 和 write2string 命令,它们用于写入而不是输出文件,可以选择一个文件或生成一个新字符串。
最后,确实重要的是不要忘记每个复合文件都应该以一个最终的 End 结束。
一旦我们运行前面的示例,我们将得到以下输出:
////////////////////////////////////////////////////////
// POTENTIAL ENERGY RESULT
////////////////////////////////////////////////////////
1 50.00 -56.486626696200
2 54.00 -56.498074637200
3 58.00 -56.505200120800
4 62.00 -56.508823168800
5 66.00 -56.509732863600 (-)
6 70.00 -56.508724734300
7 74.00 -56.506590613800
8 78.00 -56.504070086000
9 82.00 -56.501791816800
10 86.00 -56.500229017900
11 90.00 -56.499674856600 (+)
12 94.00 -56.500229018100
13 98.00 -56.501791817200
14 102.00 -56.504070082800
15 106.00 -56.506590613300
16 110.00 -56.508724733100
17 114.00 -56.509732863700 (-)
18 118.00 -56.508823172900
19 122.00 -56.505200132200
20 126.00 -56.498074642900
21 130.00 -56.486626729200
////////////////////////////////////////////////////////
随着步骤,对应步骤的角度,约束优化能量的能量加上两个极小值和一个极大值在势能中的符号。