方法文章

从ab initio热化学计算分子团簇的大气浓度

9.4K 次观看

⸱

DOI：

10.3791/60964

⸱

2020年4月8日

 ,  ,  ,  , 
1Department of Chemistry, Furman University, 2College of Charleston

本文内容

摘要

通过利用遗传算法以及半经验和ab initio量子化学的多步构型采样方法，可从低能结构的热化学性质计算弱结合分子簇的大气浓度。

摘要

大气气溶胶形成与生长的计算研究需要精确的吉布斯自由能面，该自由能面可通过气相电子结构和振动频率计算获得。这些物理量适用于其几何结构对应于势能面上极小值点的大气团簇。最小能量结构的吉布斯自由能可用于预测该团簇在不同温度和压力等条件下的大气浓度。本文介绍了一种计算成本较低的方法，该方法基于遗传算法的构型采样，随后进行一系列精度逐步提高的筛选计算。该方法首先利用半经验模型生成并演化大量构型的几何结构，然后在一系列高精度的从头算（ab initio）理论水平上对所得的独特结构进行优化。最后，对得到的最小能量结构集合计算热力学修正项，并用于计算其生成吉布斯自由能、平衡常数以及大气浓度。本文展示了该方法在环境条件下水合甘氨酸团簇研究中的应用。

引言

在气候变化的大气研究中，最不确定的参数是云粒子反射入射太阳辐射的确切程度。气溶胶是悬浮在气体中的微粒物质，可形成被称为云凝结核（cloud condensation nuclei, CCN）的云粒子，这些粒子能够散射入射辐射，从而阻止其被吸收以及由此导致的大气升温1。要深入理解这种净冷却效应，必须了解气溶胶如何增长为CCN，而这就进一步需要理解小分子团簇如何增长为气溶胶粒子。最近的研究表明，气溶胶的形成始于直径为3 nm或更小的分子团簇2；然而，利用实验技术难以探测这一尺寸范围的粒子3,4。因此，需要采用计算建模的方法来克服这一实验上的局限性。

利用下文所述的建模方法，我们可以分析任意水合团簇的生长过程。由于我们关注水在前生命环境中促进小分子形成大生物分子的作用，因此以甘氨酸为例来说明本方法。所遇到的研究挑战及所需的分析工具，与大气气溶胶和前成核团簇研究中的问题非常相似5,6,7,8,9,10,11,12,13,14,15。本文中，我们从单个孤立的甘氨酸分子出发，逐步依次添加最多五个水分子，研究其形成的水合甘氨酸团簇。最终目标是计算在室温、海平面条件下、相对湿度（RH）为100%时，大气中Gly(H2O)n=0-5团簇的平衡浓度。

这些亚纳米级分子团簇中仅有少量通过吸附其他气相分子或与已有团簇凝聚，生长为亚稳态的临界团簇（直径1-3 nm）。这些临界团簇具有有利的生长特性，可进一步形成更大的云凝结核（CCN，直径达50-100 nm），从而直接影响云的降水效率及其对入射光的反射能力。因此，深入理解分子团簇的热力学性质及其平衡分布，将有助于更准确地预测气溶胶对全球气候的影响。

气溶胶形成的描述性模型需要精确的分子团簇形成热力学数据。计算分子团簇形成的精确热力学参数需要确定最稳定的构型，这涉及寻找团簇势能面（PES）上的全局和局部极小值16。该过程称为构型采样，可通过多种技术实现，包括基于分子动力学（MD）17,18,19,20、蒙特卡洛（MC）21,22以及遗传算法（GA）23,24,25的方法。

多年来，已开发出多种协议，以在高水平理论上获得大气水合物的结构和热力学性质。这些协议在以下方面有所不同：（i）构型采样方法的选择，（ii）构型采样中所用低级别方法的性质，以及（iii）在后续步骤中用于优化结果的高级别方法的层级体系。

构型采样方法包括化学直觉26、随机采样27,28、分子动力学（MD）29,30、势阱跳跃（BH）31以及遗传算法（GA）24,25,32。与这些采样方法结合使用的最常见低层次方法是力场或半经验模型，例如 PM6、PM7 和 SCC-DFTB。随后通常采用密度泛函理论（DFT）计算，所用基组逐渐增大，并选用 Jacob 梯子更高层级上更为可靠的泛函33。在某些情况下，还会进一步采用更高层次的波函数方法，例如 MP2、CCSD(T) 以及计算成本较低的 DLPNO-CCSD(T)34,35。

Kildgaard 等人36开发了一种系统性方法，该方法在较小的水合或非水合团簇周围的斐波那契球面37上的点位添加水分子，以生成更大团簇的候选结构。基于近距离接触阈值以及不同构象体之间的均方根距离，排除不合理的和冗余的候选结构。随后采用 PM6 半经验方法进行优化，并结合密度泛函理论（DFT）和波函数方法的层级计算，以在高水平理论下获得一组低能构象体。

人工蜂群（ABC）算法38是一种新的构型采样方法，最近由Zhang等人在名为ABCluster的程序中实现，用于研究分子团簇39。Kubecka等人40使用ABCluster进行构型采样，随后采用紧束缚GFN-xTB半经验方法41进行低层级的再优化。他们进一步使用DFT方法对结构和能量进行精细优化，最终采用DLPNO-CCSD(T)方法计算能量。

无论采用何种方法，构型采样均始于在势能面上随机或非随机生成的一组分布点。每个点对应于所研究分子团簇的特定几何构型，由采样方法生成。随后通过沿梯度方向搜索，为每个点找到最近的局部极小值 "下坡" 方向上的PES。所找到的这一组极小值对应于分子簇的那些至少在一段时间内稳定的几何构型。在此，PES的形状以及表面上每一点的能量计算结果将依赖于体系的物理描述，其中更精确的物理描述会导致更高的计算代价。我们将特别采用在OGOLEM中实现的GA方法25 程序，已成功应用于多种全局优化和构型采样问题42,43,44,45，以生成初始采样点集。该PES将由PM7模型描述46 在 MOPAC2016 程序中实现47这种组合方法被采用是因为它相较于MD和MC方法能产生更多样化的点，且比对势能面（PES）进行更详细描述的方法更快地找到局部极小值。

将GA优化得到的局部极小值组合作为一系列筛选步骤的起始几何构型，从而获得一组低能量的极小值结构。该部分方案首先使用小基组的密度泛函理论（DFT）对唯一存在的GA优化结构进行优化。这组优化通常会得到更少的一组唯一局部极小值结构，相较于GA优化的半经验结构，其建模更为精确。随后，对这一较小的结构集合使用更大的基组进行新一轮DFT优化。此步骤通常会进一步减少唯一结构的数量，且相较于小基组DFT步骤，其建模更加精细。最终得到的唯一结构将被进一步优化至更严格的收敛标准，并计算其谐振振动频率。完成此步骤后，我们即获得了计算大气中团簇平衡浓度所需的全部信息。整体方法的流程在图1中以示意图形式总结。我们将采用Gaussian0949程序中实现的PW9148广义梯度近似（GGA）交换-相关泛函，以及Pople50基组的两种变体（小基组步骤使用6-31+G*，大基组步骤使用6-311++G**）。选择这一特定的交换-相关泛函与基组组合，是因其在先前研究中成功用于计算大气团簇生成吉布斯自由能的高准确性51,52。

本方案假设用户可以访问一个配备便携式批处理系统53（PBS）、MOPAC2016（http://openmopac.net/MOPAC2016.html）47、OGOLEM（https://www.ogolem.org）25、Gaussian 09（https://gaussian.com）49 和 OpenBabel54（http://openbabel.org/wiki/Main_Page）软件的高性能计算（HPC）集群，且这些软件已根据各自的安装说明完成安装。本方案中的每一步均使用一组内部开发的 shell 和 Python 2.7 脚本，这些脚本必须保存至用户 $PATH 环境变量所包含的目录中。所有运行上述程序所需的环境模块和执行权限也必须加载到用户的会话中。按照现代计算机资源标准，遗传算法代码（OGOLEM）和半经验方法代码（MOPAC）的磁盘和内存占用非常小。OGOLEM/MOPAC 的总体内存和磁盘使用量取决于所使用的线程数量，即便如此，其资源消耗仍远低于大多数 HPC 系统的处理能力。量子力学方法的资源需求则取决于团簇的大小以及所采用的理论级别。采用本方案的优势在于，用户可根据需要调整理论计算级别，以获得最终的一组低能结构，但需注意，通常计算速度越快，结果的准确性不确定性也越高。

为了表述清晰，用户的本地计算机将称为"本地计算机"，而用户可访问的高性能计算集群将称为"远程集群"。

方案

1. 寻找孤立甘氨酸和水的最低能量结构

注意：此处的目标有两个：（i）获得孤立水分子和甘氨酸分子的最低能量结构，用于遗传算法构型采样；（ii）计算这些分子在气相能量中的热力学修正值，用于大气浓度的计算。

  1. 在本地计算机上，打开 Avogadro 的新会话。
    1. 点击 Build > Insert > Peptide，在弹出的 Insert Peptide 窗口中选择 Gly，以在可视化窗口中生成一个甘氨酸单体。
    2. 点击 Extensions > Gaussian，在文本框中将第一行修改为“# pw91pw91/6-311++G** int(Acc2E=12,UltraFine) scf(conver=12) opt(tight,maxcyc=300) freq”，然后点击 Generate，并将输入文件保存为 glycine.com。
    3. 请注意，如果分子具有显著的构象柔性（如甘氨酸）55，则必须进行构象分析，以确定全局能量最低结构及其他低能构象体。OpenBabel54 提供了基于不同算法和快速力场的可靠构象搜索工具。尽管在遗传算法（GA）及后续计算过程中允许构象体松弛和相互转化，但有时仍需运行多次 GA 计算，每次从不同的初始构象体开始。
  2. 在本地计算机上，打开 Avogadro 的新会话。
    1. 点击 Build > Insert > Fragment，在 Insert Fragment 窗口中搜索“water”，以获取水分子的坐标。
    2. 点击 Extensions > Gaussian，在文本框中将第一行修改为“# pw91pw91/6-311++G** int(Acc2E=12,UltraFine) scf(conver=12) opt(tight,maxcyc=300) freq”，然后点击 Generate，并将输入文件保存为 water.com。
  3. 将两个 .com 文件传输至远程集群。登录远程集群后，通过批处理提交脚本调用 Gaussian 09 以启动计算。计算完成后，调用 OpenBabel 提取最低能量结构的笛卡尔坐标（.xyz 文件）。对于甘氨酸，执行的命令为：
    obabel -ig09 glycine.log -oxyz > glycine.xyz
    这两个 .xyz 文件将在下一步的 GA 构型采样中使用。

2. 基于遗传算法的Gly(H2O)n=1-5团簇构型采样

注意：此处的目标是使用 MOPAC47 中实现的 PM746 模型，在低成本的半经验理论水平上获得 Gly(H2O)n=1-5 的一组低能结构。工作目录必须具有与图2所示完全相同的组织和结构，以确保自定义的 shell 脚本和 Python 脚本能够正常运行。

  1. 将所有必要的脚本复制到远程集群，并将其位置添加到 $PATH
    1. 将所有脚本和模板文件放入一个文件夹（例如 scripts），然后将其复制到远程集群
    2. 确保所有脚本都具有可执行权限
    3. 在终端中输入以下命令，将脚本目录的位置添加到 $PATH 环境变量。脚本的默认位置设置为 $HOME/JoVE-demo/scripts，但用户也可以定义一个名为 $SCRIPTS_HOME 的环境变量，指向包含脚本的目录，并将 $SCRIPTS_HOME 添加到自己的路径中
      1. Bash shell：
        export SCRIPTS_HOME=/path/to/scripts
        export PATH=${SCRIPTS_HOME}:${PATH}
      2. Tcsh/Csh shell：
        setenv SCRIPTS_HOME /path/to/scripts
        setenv PATH ${SCRIPTS_HOME}:${PATH}
  2. 在远程集群上设置并运行遗传算法（GA）计算：
    1. 创建一个名为 gly-h2o-n 的目录，其中 n 为水分子的数量。
    2. 在 gly-h2o-n 目录下创建一个名为 GA 的子目录，用于运行遗传算法计算。
    3. 将 OGOLEM 输入文件（例如 pm7.ogo）、单体笛卡尔坐标文件（例如 glycine.xyz、water.xyz）以及 PBS 批处理提交脚本（例如 run.pbs）复制到 GA 目录中。
    4. 对 OGOLEM 输入文件和批处理提交文件进行必要的修改。
    5. 提交计算任务。当计算开始后，OGOLEM 将在 GA 目录中创建一个以 OGOLEM 输入文件前缀命名的新目录（例如 pm7），并将新生成的坐标文件存储在该目录中。
  3. 计算完成后，整理能量和转动常数，并利用这些信息确定哪些是唯一的低能结构：
    1. 切换到 gly-h2o-n/GA/pm7 目录并
    2. 使用以下命令提取能量并计算 GA 优化后团簇的转动常数：
      getRotConsts-GA.csh N 0 99
      其中 N 是分子团簇中的原子数，“0 99”表示 GA 种群大小为 100，索引范围从 0 到 99。该命令将生成一个名为 rotConstsData_C 的文件，其中包含所有 GA 优化团簇构型的排序列表、其能量及其转动常数。
    3. 执行以下命令：
      similarityAnalysis.py pm7 rotConstsData_C
      其中 pm7 将用作文件命名标签，用于查找并保存唯一的 GA 优化团簇。该命令将生成一个名为 uniqueStructures-pm7.data 的文件，其中包含唯一 GA 优化构型的排序列表。这是 Gly(H2O)n 团簇在 PM7 理论水平上优化得到的唯一局部极小结构列表，这些结构现已准备好使用 DFT 进一步优化。
  4. 返回到 gly-h2o-n/GA 目录，并使用 combine-GA.csh 脚本合并多个可比的 GA 运行结果。语法为：
    combine-GA.csh

3. 使用小基组的量子力学方法进行优化

注意：此处的目标是利用更精确的量子力学描述方法，优化 Gly(H2O)n=1-5 团簇的构型采样，从而获得一组更小但更准确的 Gly(H2O)n=1-5 团簇结构。此步骤的初始结构来自第 2 步的输出结果。

  1. 准备并运行小基组DFT计算：
    1. 在gly-h2o-n目录下创建一个名为QM的子目录。在QM目录下，再创建一个名为pw91-sb的子目录。
    2. 将唯一结构列表（uniqueStructures-pm7.data）从gly-h2o-n/GA目录复制到QM/pw91-sb目录中。
    3. 切换当前目录至gly-h2o-n/QM/pw91-sb。
    4. 使用以下命令运行小基组DFT构型采样脚本：
      run-pw91-sb.csh uniqueStructures-pm7.data sb QUEUE 10
      其中sb是该组计算的标签，QUEUE是计算集群上首选的队列名称，10表示将10个计算任务合并为一个批处理作业。该脚本将自动为Gaussian 09生成输入文件并提交所有计算任务。若输入“test”作为“QUEUE”参数，则执行一次试运行。
  2. 当提交的计算任务完成后，提取并分析结果。
    1. 使用以下命令提取能量，并计算经小基组优化后的团簇的转动常数：
      getRotConsts-dft-sb.csh pw91 N
      其中pw91表示使用了PW91密度泛函，N为团簇中的原子数。该命令将生成一个名为rotConstsData_C的文件。
    2. 接下来使用以下命令识别唯一结构：
      similarityAnalysis.py sb rotConstsData_C
      其中sb用作文件命名标签。此时，将在文件uniqueStructures-sb.data中保存一份在PW91/6-31+G*理论水平下优化得到的唯一构型列表。
  3. 返回至gly-h2o-n/QM目录，并使用combine-QM.csh脚本合并多个可比的量子力学（QM）计算结果。其语法为：
    combine-QM.csh

4. 使用大基组的量子力学方法进行进一步优化

注意：此处的目标是利用更精确的量子力学描述，进一步优化 Gly(H2O)n=1-5 团簇的构型采样。此步骤的起始结构来自第 3 步的输出结果。

  1. 使用更大的基组提交更可靠的计算。
    1. 在 QM 目录下创建一个名为 pw91-lb 的子目录。
    2. 将唯一结构列表（uniqueStructures-sb.data）从 gly-h2o-n/QM 目录复制到 gly-h2o-n/QM/pw91-lb 目录，并切换至该目录。
    3. 运行大基组DFT构型采样脚本，命令为：
      run-pw91-lb.csh uniqueStructures-sb.data lb QUEUE 10
      其中 lb 是此组计算的标签，QUEUE 是计算集群上首选的队列，10 表示每批作业中包含 10 次计算。该脚本将自动为 Gaussian 09 生成输入文件并提交所有计算任务。若将 'QUEUE' 设为 'test'，则进行一次试运行测试。
  2. 当提交的计算完成后，提取并分析数据
    1. 使用以下命令计算大基组优化团簇的转动常数：
      getRotConsts-dft-lb.csh pw91 N
      其中 pw91 表示使用了 PW91 密度泛函，N 表示团簇中的原子数。
    2. 现在使用以下命令识别唯一结构：
      similarityAnalysis.py lb rotConstsData_C
      其中 lb 用作文件命名标签。此时您已获得一份在 PW91/6-311++G** 理论水平上优化后的唯一构型列表，保存在文件 uniqueStructures-lb.data 中。

5. 最终能量与热力学校正计算

注意：此处的目标是使用大型基组和超精细积分网格，获得Gly(H2O)n=1-5团簇的振动结构和能量，以计算所需的热化学校正。

  1. 基于前一步的结果，提交更可靠的计算。
    1. 在 QM/pw91-lb 目录下创建一个名为 ultrafine 的子目录。然后将唯一结构列表（uniqueStructures-lb.data）从 QM/pw91-lb 目录复制到 QM/pw91-lb/ultrafine 目录，并切换至该目录。
    2. 使用以下命令提交超细网格大基组 DFT 计算脚本：
      run-pw91-lb-ultrafine.csh uniqueStructures-lb.data uf QUEUE 10
      其中 uf 是本组计算的标签，QUEUE 是计算集群上首选的队列，10 表示每批作业中包含 10 个计算任务。该脚本将自动为 Gaussian 09 生成输入文件并提交所有计算任务。若将 'QUEUE' 替换为 'test'，则执行一次试运行测试。
  2. 当提交的计算完成后，提取并分析数据
    1. 使用以下命令提取能量，并计算大基组优化团簇的转动常数：
      getRotConsts-dft-lb-ultrafine.csh pw91 N
      其中 pw91 表示使用了 PW91 密度泛函，N 是团簇中的原子数目。
    2. 现在使用以下命令识别唯一结构：
      similarityAnalysis.py uf rotConstsData_C
      其中 uf 用作文件命名标签。此时您将获得一份在 PW91/6-311++G** 理论水平上优化后的唯一构型列表，保存在文件 uniqueStructures-uf.data 中。
  3. 最终提取用于计算热力学修正所需的信息，并利用该信息计算热力学修正值。
    1. 提取最终的电子能量、转动常数和振动频率，并使用以下命令计算热力学修正值：
      run-thermo-pw91.csh uniqueStructures-uf.data
    2. 将命令行输出内容复制粘贴到名为 'gly-h2o-n.xlsx' 的 Excel 工作簿中的 'Raw_Energies' 工作表。您需要对单体（甘氨酸和水）以及每种水合物中能量最低的构型（gly-h2o-n，其中 n=1,2,…）执行相同操作。
    3. 当原始能量被添加至 'gly-h2o-n.xlsx' 工作簿的第一个工作表后，后续的 'Binding_Energies' 和 'Hydrate_Distribution' 工作表将自动更新。特别是 'Hydrate_Distribution' 工作表会给出在不同温度（例如 298.15 K）、相对湿度（20%、50%、100%）以及水（[H2O]）和甘氨酸（[Glycine]）初始浓度条件下，各水合物的平衡浓度。这些计算背后的理论将在下一步中描述。

6. 计算室温下海平面处Gly(H2O)n=0-5团簇的大气浓度

注意：首先将上一步生成的热力学数据复制到电子表格中，计算逐级水合的吉布斯自由能。然后利用吉布斯自由能计算每一步水合反应的平衡常数。最后，求解一组线性方程，得到在给定单体浓度、温度和压力条件下的各种水合物的平衡浓度。

  1. 首先建立甘氨酸逐级水合反应的化学平衡体系，如下所示：
    Hydration reactions of Gly in water; chemical equations; equilibrium process.
  2. 计算平衡常数 Kn 使用 Kn = e-ΔGn/(kBT)，其中 n 是水合程度，ΔGn 是反应步骤 n 的吉布斯自由能变化th 水合反应， kB 是玻尔兹曼常数，以及 T 是温度。
    Equilibrium equations, glycol and water interaction, chemical formula, scientific analysis
  3. 建立质量守恒方程，假设水合与非水合甘氨酸团簇的平衡浓度之和等于孤立甘氨酸的初始浓度 [Gly]0将这组包含六个联立方程的系统通过平衡常数表达式的代数重排，改写为
    Equilibrium equations with glycine-water complexes, depicting molecular binding states.
  4. 解上述方程组，得到甘氨酸（Gly(H）的平衡浓度2On = 0-5 使用实验值56,57,58 大气中甘氨酸的浓度，[Gly]0 = 2.9 × 10⁶ cm-3，以及在相对湿度为100%、温度为298.15 K时大气中水的浓度59，[H2O] = 7.7 × 1017 厘米-3.

结果

本实验方案得到的第一组结果应为通过构型采样过程获得的Gly(H2O)n=1-5的一系列低能结构。这些结构已在PW91/6-311++G**理论水平上完成优化，并被认为在本文研究目的下具有足够的准确性。目前没有证据表明PW91/6-311++G**会系统性地低估或高估这些团簇的结合能。该方法在预测结合能方面相对于MP2/CBS32和[DLPNO-]CCSD(T)/CBS60,61的估计值以及实验值52表现出较大的波动性，大多数其他密度泛函也存在类似情况。通常情况下，每个n = 1 – 5的体系都应产生若干个能量接近最低能量结构约5 kcal mol-1范围内的低能构型。为简洁起见，本文仅关注由run-thermo-pw91.csh脚本生成的第一个结构。图3展示了Gly(H2O)n=0-5团簇中电子能量最低的异构体。可以看出，随着水分子数量的增加，氢键网络的复杂性逐渐增强，甚至在n = 5时从主要呈平面的结构演变为三维笼状结构。本文后续部分所使用的能量和热力学量均对应于这五个特定团簇。

表1 包含执行本实验方案所需的热力学参数。表2 展示了 run-thermo-pw91.csh 脚本输出结果的一个示例，其中列出了电子能、振动零点能校正项以及在三个不同温度下的热力学校正项。对于每个团簇（每一行），E[PW91/6-311++G**] 表示在PW91/6-311++G**理论水平下、使用超精细积分网格计算得到的气相电子能，单位为哈特里（Hartree），以及零点振动能（ZPVE），单位为kcal mol-1。在每个温度（216.65 K、273.15 K 和 298.15 K）下，列出了相应的热力学校正项：∆H 为生成焓，单位为kcal mol-1；S 为生成熵，单位为cal mol-1；∆G 为吉布斯自由能变，单位为kcal mol-1。表3 展示了总水合吉布斯自由能变以及逐级水合过程的计算示例。以下为反应的总水合吉布斯自由能变的一个计算示例：

figure-results-1

从计算电子能量 EPW91 开始

figure-results-2

其中 EPW91[糖原∙(H2O)] 来自于 表2 C列和E列PW91[Gly] 和 EPW91[H2O] 来自于 表1 B列。接下来我们计算气相总能量变化 ΔE(0) 通过包含反应零点振动能的变化

figure-results-3

以获得D列数据。其中，ΔEPW91/6−311++G**取自表3的C列，EZPVE[Gly ∙ (H2O)]取自表2的D列，EZPVE[Gly]和EZPVE[H2O]取自表1的C列。为简洁起见，我们直接进入室温下的团簇分析，因此跳过216.65 K和273.15 K的数据。在室温条件下，我们通过校正气相能量变化来计算反应的焓变ΔH，其表达式为

figure-results-4

其中 ΔE(0) 来自于 表3 D 列， ΔH[糖原∙(H2O)] 来自于 表2 第 K 列，以及 ΔH[甘氨酸]和 ΔH[H2O] 来自 表1 第 J 列。最后，我们计算该反应的吉布斯自由能变化 ΔG 作为

figure-results-5

其中 ΔH 取自 表3 第 I 列 S[糖原∙(H2O)] 来自于 表2 第 L 列，以及 S[甘氨酸]和 S[H2O] 来自于 表1 第K列。请注意，此处的熵值必须转换为 kcal/mol 单位-1 K-1 在此步骤中。

我们现在已具备必要的数据量，可按照步骤6所示计算水合甘氨酸的大气浓度。结果应与表4中所示数据相似，但可能存在微小的数值差异。表4展示了通过将步骤6.2中的六个方程组整合为一个矩阵方程并求解后得到的平衡态水合物浓度。我们首先确认该方程组可表示为

figure-results-6

其中 Kn 是n的平衡常数th 甘氨酸的逐步水合 w 是大气中水的浓度， g 是大气中分离出的甘氨酸的初始浓度，以及 gn 是Gly(H)的平衡浓度2On如果我们将上述方程重写为 Ax = b，我们得到 x = A−1b 其中 A−1 矩阵的逆 A。该逆矩阵可使用电子表格内置函数轻松计算，如下所示 表4 以获得最终结果。

图4展示了在100%相对湿度和1个大气压下，水合甘氨酸的平衡浓度随温度变化的情况，该数据基于表4的计算结果。结果显示，当温度从298.15 K降低至216.65 K时，无水甘氨酸（n=0）的浓度减小，而水合甘氨酸的浓度则增加。特别是甘氨酸二水合物（n=2）的浓度随温度降低显著上升，而其他水合物浓度的变化则相对较小。这种温度与水合物浓度之间的负相关关系符合预期，即在较低温度下水合过程的吉布斯自由能更低，有利于水合物的形成。

图5展示了在298.15 K和1个大气压下，甘氨酸水合物平衡浓度对相对湿度的依赖性。图中清楚表明，当相对湿度从20%升高至100%时，水合物（n>0）的浓度增加，而未水合甘氨酸（n=0）的浓度相应减少。相对湿度与水合物浓度之间的直接相关性再次表明，在较高相对湿度下，更多的水分子有利于水合物的形成。

如上所述，本方案可对大气中水合甘氨酸的种群分布提供定性理解。假设孤立甘氨酸的初始浓度为每立方厘米290万个分子，则除T=216.65 K且相对湿度（RH）为100%的条件外，在大多数条件下无水甘氨酸（n=0）均为最丰富的物种。二水合物（n=2）在三个温度下均具有最低的逐级吉布斯自由能水合值，是在本研究考虑的条件下最丰富的水合物。一水合物（n=1）和更大尺寸的水合物（n≥3）预计含量可忽略不计。通过观察图3，n = 1–4团簇的丰度可与团簇中氢键网络的稳定性和应变相关联。这些团簇中的水分子以接近多种氢键环状结构的几何构型与甘氨酸的羧酸基团形成氢键，从而表现出特别高的稳定性。

figure-results-7
图 1：当前流程的示意图。 由遗传算法（GA）生成的大量候选结构池，经过一系列 PW91 几何优化后逐步优化，直至获得一组收敛的结构。计算这些结构的振动频率，并用于计算其生成吉布斯自由能，进而用于计算在环境条件下团簇的平衡浓度。请点击此处查看此图的放大版本。

figure-results-8
图2：每个簇的代表性目录结构。 本方案中包含的内部脚本要求采用如上所示的目录结构，其中 n 为水分子的数量。对于 gly-h2o-n 中的每一个 n，均包含以下子目录：GA 用于遗传算法，内含 GA/pm7 目录；QM 用于量子力学计算，包含 QM/pw91-sb 目录（对应 PW91/6-31+G* 方法）、QM/pw91-lb 目录（对应 PW91/6-311++G** 方法），以及 QM/pw91-lb/ultrafine 目录，用于在超精细积分网格上进行结构优化和最终的振动频率计算。请点击此处查看该图的放大版本。

figure-results-9
图3：Gly(H2O)n=0-5 的代表性低能结构。 这些团簇是在PW91/6-311++G**理论水平上优化得到的电子能全局最小值结构。请点击此处查看该图的放大版本。

figure-results-10
图 4：Gly(H2O)n=0-5 在 100% 相对湿度和 1 atm 压力下的温度依赖性。 水合物的浓度以分子数每立方厘米（molecules cm-3）为单位给出。请点击此处查看此图的放大版本。

figure-results-11
图5：298.15 K 和 1 atm 压力下 Gly(H2O)n=0-5 的相对湿度依赖性。 水合物的浓度以分子数每立方厘米（molecules cm-3）为单位给出。请点击此处查看此图的放大版本。

E[PW91/6-311++G**]216.65 K273.15 K298.15 K
LB-UFZPVE∆HS∆G∆HS∆G∆HS∆G
水-76.43050013.041.7242.595.542.1744.443.082.3745.141.96
甘氨酸-284.43483848.552.6569.5336.143.7073.8132.094.2275.6130.22

表1：单体能量。 电子能量的单位为哈特里（Hartree），其他所有量的单位均为千卡每摩尔（kcal mol-1）。水和甘氨酸在PW91/6-311++G**理论水平上进行了优化，并计算了振动频率。在1 atm压力和298.15 K温度下的热力学修正值使用thermo.pl脚本计算得到。

E[PW91/6-311++G**]0 K216.65 K273.15 K298.15 K
n名称LB-UF零点振动能 (ZPVE)∆HS∆G∆HS∆G∆HS∆G
1gly-h2o-1-360.8848163.963.6180.1250.225.1286.2745.525.8588.8343.33
2gly-h2o-2-437.3376379.334.5390.8664.176.4698.7858.817.40102.0656.30
3gly-h2o-3-513.7862094.525.67105.0877.428.08114.9471.199.23119.0068.27
4gly-h2o-4-590.23667109.806.03104.9891.308.78116.2184.4010.11120.8781.14
5gly-h2o-5-666.68845125.807.26121.70106.6910.47134.8399.4412.01140.2496.00

表2：簇能量 最低能量Gly(H2On=1-5 使用我们所述程序发现的结构 图1电子能量的单位为哈特里（Hartree），其他所有物理量的单位均为千卡每摩尔（kcal mol）-1.

总水合：Gly + nH2O <-> Gly(H2O)n逐级水合：Gly(H2O)n-1 + H2O <-> Gly(H2O)n
E[PW91/6-311++G**]216.65273.15298.15216.65273.15298.15
n体系名称LB-UF∆E(0)∆H(T)∆G(T)∆H(T)∆G(T)∆H(T)∆G(T)LB-UF∆E(0)∆H(T)∆G(T)H(T)∆G(T)∆H(T)∆G(T)
1gly-h2o-1-12.22-9.85-10.61-3.68-10.61-1.87-10.59-1.07-12.22-9.85-10.61-3.68-10.61-1.87-10.59-1.07
2gly-h2o-2-26.22-21.53-23.10-9.27-23.11-5.66-23.09-4.06-14.00-11.68-12.49-5.59-12.50-3.79-12.50-2.99
3gly-h2o-3-37.56-30.72-32.88-12.90-32.87-7.69-32.82-5.38-11.34-9.19-9.78-3.63-9.76-2.03-9.73-1.32
4gly-h2o-4-50.10-40.34-43.48-15.87-43.54-8.71-43.51-5.55-12.54-9.62-10.60-2.97-10.67-1.02-10.69-0.17
5gly-h2o-5-63.45-51.41-55.42-20.58-55.51-11.48-55.48-7.45-13.35-11.07-11.94-4.71-11.97-2.77-11.97-1.90

表3：水合能。 Gly(H2O)n=1-5 的总水合能及逐级水合能在 kcal mol-1 单位下的数值。其中，E[PW91/6-311++G**] 表示电子能的变化，∆E(0) 表示经零点振动能（ZPVE）校正后的能量变化，∆H(T) 表示在温度 T 下的焓变，∆G(T) 表示每个 Gly(H2O)n=1-5 团簇的水合吉布斯自由能变化。

平衡水合物分布随温度和相对湿度的变化 
T=298.15KT=273.15KT=216.65K
Gly(H2O)nRH=100%RH=50%RH=20%RH=100%RH=50%RH=20%RH=100%RH=50%RH=20%
01.3E+062.2E+062.7E+061.1E+062.0E+062.7E+066.1E+051.5E+062.5E+06
12.3E+051.9E+059.5E+042.0E+051.9E+059.9E+041.2E+051.5E+059.5E+04
21.0E+064.3E+058.4E+041.3E+066.1E+051.3E+051.8E+061.1E+063.0E+05
32.8E+055.8E+044.5E+033.2E+057.4E+046.3E+033.1E+059.6E+041.0E+04
41.1E+041.1E+033.4E+011.3E+041.5E+035.0E+011.1E+041.8E+037.5E+01
57.5E+033.9E+024.9E+001.2E+047.2E+029.7E+002.4E+041.9E+033.1E+01

表4：Gly(H2O)n=0-5 在不同温度（T=298.15 K、273.15 K、216.65 K）和相对湿度（RH=100%、50%、20%）下的平衡水合物浓度。 水合物浓度的单位为分子数每立方厘米（molecules cm-3），计算基于实验值56,57,58,，其中[Gly]0 = 2.9 × 106 cm-3，而[H2O]在相对湿度为100%时，在T = 298.15 K、273.15 K和216.65 K下分别为7.7 × 1017 cm-3、1.6 × 1017 cm-3和9.9 × 1014 cm-359。

补充文件。 请点击此处下载这些文件。

讨论

本方案生成数据的准确性主要取决于三个方面：（i）步骤2所采样的构型多样性，（ii）体系电子结构的准确性，（iii）热力学修正的准确性。上述每个因素均可通过修改所附脚本中的方法加以改进。第一个因素可通过使用更大规模的随机生成结构初始池、增加遗传算法（GA）的迭代次数以及放宽GA中相关判据的定义来有效克服。此外，也可采用其他半经验方法，例如自洽电荷密度泛函紧束缚（SCC-DFTB）62模型和有效片段势（EFP）63模型，以探究不同物理描述带来的影响。该方法的主要局限在于无法形成或断裂共价键，这意味着单体结构被冻结。遗传算法过程仅根据半经验描述，寻找这些被冻结单体之间最稳定的相对位置。

可以通过多种方式提高体系电子结构的准确性，每种方法都有其相应的计算成本。可以选择更优的密度泛函，例如 M06-2X64 和 wB97X-V65，或更精确的量子力学（QM）方法，例如 Møller-Plesset66,67,68（MPn）微扰理论和耦合簇69（CC）方法，以改进对体系的物理描述。在泛函的层级结构中，性能通常随着从广义梯度近似（GGA）泛函（如 PW91）向长程校正杂化泛函（如 wB97X-D）和 meta-GGA 杂化泛函（如 M06-2X）过渡而逐步提升。

DFT 方法的缺点是无法系统地收敛到精确值；然而，DFT 方法在计算上成本较低，并且存在多种适用于不同应用的功能泛函。

采用波函数方法（如MP2和CCSD(T)）结合相关一致基组（其基组基数依次增加，例如[aug-]cc-pV[D,T,Q,...]Z）计算得到的能量会系统性地收敛至完备基组极限，但随着体系尺寸增大，每次计算的计算成本变得难以承受。通过使用显式关联基组70以及外推至完备基组（CBS）71极限，可进一步优化电子结构。我们近期的研究表明，采用密度拟合的显式关联二阶Møller-Plesset微扰法（DF-MP2-F12）所获得的能量接近MP2/CBS计算的结果32。若要修改当前协议以采用不同的电子结构方法，需进行两个步骤：(i) 根据软件规定的语法准备模板输入文件；(ii) 编辑 run-pw91-sb.csh、run-pw91-lb.csh 和 run-pw91-lb-ultrafine.csh 脚本，以生成符合该软件要求的正确输入文件语法及正确的提交脚本。

最后，热力学修正的准确性取决于电子结构方法以及对全局极小值附近势能面（PES）的描述。准确描述势能面需要计算核自由度位移相对于势能面的三阶及更高阶导数，例如四次力场72,73（QFF），这是一项计算成本极高的任务。当前方案采用谐振子近似来处理振动频率，因此只需计算势能面的二阶导数。然而，对于高度非谐性的体系（如非常柔性的分子和对称双势阱），由于真实势能面与谐振势能面之间存在显著差异，该方法会遇到问题。此外，若采用计算代价较高的高精度电子结构方法获得高质量势能面，将进一步加剧振动频率计算的成本问题。一种解决方法是结合高质量电子结构计算得到的电子能与在较低精度势能面上计算的振动频率，从而在计算成本与准确性之间取得平衡。当前方案可按照前文所述方式修改以采用不同的势能面描述；此外，也可直接修改脚本和模板中的振动频率关键词，以计算非谐性振动频率。

对于任何构型采样方案而言，两个关键问题分别是势能面的初始采样方法以及用于识别每个团簇的判据。在我们以往的研究中已广泛采用了多种方法。针对第一个问题，即势能面的初始采样方法，我们基于以下因素选择了使用遗传算法（GA）结合半经验方法。化学直觉引导的构型采样26、随机采样以及分子动力学（MD）29,30在处理超过10个单体的团簇时，通常难以稳定地找到可能的全局极小结构，这一点在我们对水团簇的研究中已观察到18。我们曾成功应用盆地跳跃法（BH）研究（H2O）11复杂的势能面74，但该方法需要手动引入一些BH算法未能发现的潜在低能异构体。BH与GA在寻找水团簇（H2O）n=10-20全局极小结构方面的性能比较表明，GA始终比BH更快地找到全局极小值75。OGOLEM和CLUSTER中实现的GA具有高度通用性，因为它可应用于任意分子团簇，并能与大量支持经典力场、半经验方法、密度泛函以及从头算（ab initio）的计算程序对接。选择PM7方法是因其计算速度快且具备合理的准确性；实际上，任何其他半经验方法都会带来显著更高的计算成本。

关于第二个问题，我们探索了使用不同的标准来识别独特结构，这些标准包括电子能量、偶极矩、重叠均方根偏差（RMSD）和转动常数。使用偶极矩存在困难，因为偶极矩的各个分量依赖于分子的取向，而总偶极矩对几何结构的差异极为敏感，导致难以设定判断结构是否相同或独特的阈值。实践证明，结合电子能量和转动常数的方法最为有效。

目前判断两种结构是否唯一的标准基于能量差阈值0.10 kcal mol-1和转动常数差异1%。因此，当两种结构的能量差超过0.10 kcal mol-1（约0.00015 a.u.）且其三个转动常数（A、B、C）中任意一个的差异超过1%时，即视为不同结构。多年来的大量内部基准测试表明，这些阈值是合理的选择。我们的构型采样方法和筛选策略已应用于多种体系，包括弱结合的多环芳烃-水复合物76,77，以及强结合的含氨和胺类的三元硫酸盐水合物32。对于存在多种质子化状态的团簇，最佳策略是进行多个遗传算法（GA）计算，每次计算从不同质子化状态的单体开始，以确保充分考虑不同质子化状态的结构。然而，低层级的DFT计算通常允许在几何优化过程中发生质子化状态的转变，从而无论初始几何构型如何，最终均会得到最稳定的质子化状态。

我们的遗传算法（GA）构型采样方法即使对于柔性分子也应能良好适用，前提是GA程序与通用的、非参数化的方法相耦合，从而允许单体在GA运行过程中采取不同的构型。例如，将GA与PM7方法耦合可允许单体结构发生变化，但如果其化学键断裂（如质子化状态改变时可能发生的情况），这些结构可能会因不符合要求而被舍弃。

我们已考虑了多种修正谐波近似缺陷的方法，尤其是由低振动频率引起的那些问题。将准谐波近似引入当前方法并不困难。然而，准谐波方法本身仍存在一些疑问，特别是关于应用该方法时所采用的截断频率的设定问题。此外，尽管传统观点认为准RRHO近似应优于RRHO近似，但目前仍缺乏严格的基准研究来检验准RRHO近似的可靠性。

因此，本方案可推广至任何非共价结合的气相分子团簇体系。通过修改脚本和模板，也可推广至使用任何半经验方法、电子结构计算方法及软件、以及振动分析方法及软件。这要求使用者熟悉 Linux 命令行界面、Python 脚本编写和高性能计算。Linux 操作系统的陌生语法和界面，以及缺乏脚本编写经验，是本方案中最大的障碍，也是新学生最易遇到困难之处。该方案已在本课题组多年来的多种应用中成功使用，主要聚焦于硫酸和氨对气溶胶形成的影响。对该方案的进一步改进将包括与更多电子结构计算软件的更稳健接口、遗传算法的替代实现方式，以及可能采用更新的方法以更快速地计算电子和振动能。我们目前的应用正在探究氨基酸在当前大气中气溶胶形成的早期阶段，以及在前生命环境中形成较大生物分子过程中的重要性。

致谢

本项目得到了美国国家科学基金会（GCS）的资助，项目编号为 CHE-1229354、CHE-1662030、CHE-1721511 和 CHE-1903871，以及阿诺德和梅布尔·贝克曼基金会贝克曼学者奖（AGG）和巴里·M·戈德华特奖学金（AGG）的支持。研究使用了水星联盟（MERCURY Consortium）的高性能计算资源（http://www.mercuryconsortium.org）。

材料

本文使用的材料清单
姓名公司目录编号评论
Avogadrohttps://avogadro.cc开源分子可视化程序
Gaussian [09/16] 软件http://www.gaussian.com/商业从头算电子结构程序
MOPAC 2016http://openmopac.net/MOPAC2016.html开源半经验程序
OGOLEM 软件https://www.ogolem.org基于遗传算法的全局优化程序
OpenBabelhttp://openbabel.org/wiki/Main_Page开源化学信息学库
calcRotConsts.pyShields课题组，化学系，弗曼大学计算转动常数的 Python 脚本
calcSymmetry.cshShields课题组，化学系，弗曼大学用于根据笛卡尔坐标计算分子对称数的 Shell 脚本
combine-GA.cshShields课题组，化学系，弗曼大学用于合并来自不同 GA 目录的能量和转动常数的 Shell 脚本
combine-QM.cshShields课题组，化学系，弗曼大学用于合并来自不同量子力学目录的能量和转动常数的 Shell 脚本
gaussianE.cshShields课题组，化学系，弗曼大学用于提取 Gaussian 09 能量的 Shell 脚本
gaussianFreqs.cshShields课题组，化学系，弗曼大学用于提取 Gaussian 09 振动频率的 Shell 脚本
getrotconstsShields课题组，化学系，福尔曼大学根据分子的笛卡尔坐标计算转动常数的可执行程序
getRotConsts-dft-lb.cshShields课题组，化学系，弗曼大学用于计算大批量大基组DFT优化结构的转动常数的Shell脚本
getRotConsts-dft-lb-ultrafine.cshShields课题组，化学系，福尔曼大学用于计算一批超精细DFT优化结构的转动常数的Shell脚本
getRotConsts-dft-sb.cshShields课题组，化学系，弗曼大学用于计算一批小基组DFT优化结构的转动常数的Shell脚本
getRotConsts-GA.cshShields课题组，化学系，弗曼大学用于计算一批遗传算法优化结构的转动常数的 Shell 脚本
global-minimum-coords.xyzShields课题组，化学系，福尔曼大学糖-(h2o)n全局最小结构的笛卡尔坐标，其中n=0-5
make-thermo-gaussian.cshShields课题组，化学系，福尔曼大学从高斯输出文件中提取数据并为 thermo.pl 脚本生成输入文件的 Shell 脚本
ogolem-input-file.ogoShields课题组，化学系，福尔曼大学Ogolem 样本输入文件
ogolem-submit-script.pbsShields课题组，化学系，福尔曼大学用于Ogolem计算的PBS批量提交文件
README.docxShields课题组，化学系，弗曼大学帮助读者有效使用脚本的说明
runogolem.cshShields课题组，化学系，福尔曼大学运行 OGOLEM 的 Shell 脚本
run-pw91-lb.cshShields课题组，化学系，福尔曼大学用于运行大批量大基组DFT优化计算的Shell脚本
run-pw91-lb-ultrafine.cshShields课题组，化学系，弗曼大学用于运行一批超精细DFT优化计算的Shell脚本
run-pw91-sb.cshShields课题组，化学系，弗曼大学用于运行小基组DFT优化计算批处理的Shell脚本
run-thermo-pw91.cshShields课题组，化学系，福尔曼大学用于计算一批 DFT 优化结构的热力学修正的 Shell 脚本
similarityAnalysis.pyShields课题组，化学系，弗曼大学基于转动常数和能量确定唯一结构的Python脚本
对称性Shields课题组，化学系，福尔曼大学根据笛卡尔坐标计算分子对称性的可执行程序
对称性.c(C) 1996, 2003 S. Patchkovskii, Serguei.Patchkovskii@sympatico.ca用于根据笛卡尔坐标确定分子对称性的 C 代码
template-marcy.pbsShields课题组，化学系，弗曼大学使用 OGOLEM 的 PBS 提交脚本模板
template-pw91.comShields课题组，化学系，福尔曼大学模板 Gaussian 09 输入文件
template-pw91-HL.comShields课题组，化学系，福尔曼大学超精细DFT优化的Gaussian 09输入模板
thermo.plhttps://www.nist.gov/mml/csd/化学信息学研究组/产品与服务/理想气体计算程序用于计算理想气体热力学校正的 Perl 开源脚本
gly-h2o-n.xlsxShields课题组，化学系，弗曼大学完整实验方案的 Excel 表格
table-1.xlsxShields课题组，化学系，福尔曼大学Excel 电子表格
table-2.xlsxShields课题组，化学系，弗曼大学Excel 电子表格
table-3.xlsxShields课题组，化学系，弗曼大学Excel 电子表格
table-4.xlsxShields课题组，化学系，弗曼大学Excel 电子表格
水.xyzShields课题组，化学系，弗曼大学水的笛卡尔坐标
glycine.xyzShields课题组，化学系，弗曼大学甘氨酸的笛卡尔坐标

参考文献

  1. Foster, P., Ramaswamy, V. Climate Change 2007 The Scientific Basis. Solomon, S., Qin, D., Manning, M., Chen, Z., Marquis, M., Averyt, K. B., Tignor, M., Miller, H. L. , Cambridge University Press. Cambridge, U.K. (2007).
  2. Kulmala, M., et al. Toward direct measurement of atmospheric nucleation. Science. 318 (5847), 89-92 (2007).
  3. Sipila, M., et al. The role of sulfuric acid in atmospheric nucleation. Science. 327 (5970), 1243-1246 (2010).
  4. Jiang, J., et al. First measurement of neutral atmospheric cluster and 1 - 2 nm particle number size distributions during nucleation events. Aerosol Science and Technology. 45 (4), (2011).
  5. Dunn, M. E., Pokon, E. K., Shields, G. C. Thermodynamics of forming water clusters at various Temperatures and Pressures by Gaussian-2, Gaussian-3, Complete Basis Set-QB3, and Complete Basis Set-APNO model chemistries; implications for atmospheric chemistry. Journal of the American Chemical Society. 126 (8), 2647-2653 (2004).
  6. Pickard, F. C., Pokon, E. K., Liptak, M. D., Shields, G. C. Comparison of CBSQB3, CBSAPNO, G2, and G3 thermochemical predictions with experiment for formation of ionic clusters of hydronium and hydroxide ions complexed with water. Journal of Chemical Physics. 122, 024302(2005).
  7. Pickard, F. C., Dunn, M. E., Shields, G. C. Comparison of Model Chemistry and Density Functional Theory Thermochemical Predictions with Experiment for Formation of Ionic Clusters of the Ammonium Cation Complexed with Water and Ammonia; Atmospheric Implications. Journal of Physical Chemistry A. 109 (22), 4905-4910 (2005).
  8. Alongi, K. S., Dibble, T. S., Shields, G. C., Kirschner, K. N. Exploration of the Potential Energy Surfaces, Prediction of Atmospheric Concentrations, and Vibrational Spectra of the HO2•••(H2O)n (n=1-2) Hydrogen Bonded Complexes. Journal of Physical Chemistry A. 110 (10), 3686-3691 (2006).
  9. Allodi, M. A., Dunn, M. E., Livada, J., Kirschner, K. N. Do Hydroxyl Radical-Water Clusters, OH(H2O)n, n=1-5, Exist in the Atmosphere. Journal of Physical Chemistry A. 110 (49), 13283-13289 (2006).
  10. Kirschner, K. N., Hartt, G. M., Evans, T. M., Shields, G. C. In Search of CS2(H2O)n=1-4 Clusters. Journal of Chemical Physics. 126, 154320(2007).
  11. Hartt, G. M., Kirschner, K. N., Shields, G. C. Hydration of OCS with One to Four Water Molecules in Atmospheric and Laboratory Conditions. Journal of Physical Chemistry A. 112 (19), 4490-4495 (2008).
  12. Morrell, T. E., Shields, G. C. Atmospheric Implications for Formation of Clusters of Ammonium and 110 Water Molecules. Journal of Physical Chemistry A. 114 (12), 4266-4271 (2010).
  13. Temelso, B., et al. Quantum Mechanical Study of Sulfuric Acid Hydration: Atmospheric Implications. Journal of Physical Chemistry A. 116 (9), 2209(2012).
  14. Husar, D. E., Temelso, B., Ashworth, A. L., Shields, G. C. Hydration of the Bisulfate Ion: Atmospheric Implications. Journal of Physical Chemistry A. 116 (21), 5151-5163 (2012).
  15. Bustos, D. J., Temelso, B., Shields, G. C. Hydration of the Sulfuric Acid – Methylamine Complex and Implications for Aerosol Formation. Journal of Physical Chemistry A. 118 (35), 7430-7441 (2014).
  16. Wales, D. J., Scheraga, H. A. Global optimization of clusters, crystals, and biomolecules. Science. 27 (5432), 1368-1372 (1999).
  17. Day, M. B., Kirschner, K. N., Shields, G. C. Global search for minimum energy (H2O)n clusters, n = 3 - 5. The Journal of Physical Chemistry A. 109 (30), 6773-6778 (2005).
  18. Shields, R. M., Temelso, B., Archer, K. A., Morrell, T. E., Shields, G. C. Accurate predictions of water cluster formation, (H2O)n=2-10. The Journal of Physical Chemistry A. 114 (43), 11725-11737 (2010).
  19. Temelso, B., Archer, K. A., Shields, G. C. Benchmark structures and binding energies of small water clusters with anharmonicity corrections. The Journal of Physical Chemistry A. 115 (43), 12034-12046 (2011).
  20. Temelso, B., Shields, G. C. The role of anharmonicity in hydrogen-bonded systems: The case of water clusters. The Journal of Chemical Theory and Computation. 7 (9), 2804-2817 (2011).
  21. Von Freyberg, B., Braun, W. Efficient search for all low energy conformations of polypeptides by Monte Carlo methods. The Journal of Computational Chemistry. 12 (9), 1065-1076 (1991).
  22. Rakshit, A., Yamaguchi, T., Asada, T., Bandyopadhyay, P. Understanding the structure and hydrogen bonding network of (H2O)32 and (H2O)33: An improved Monte Carlo temperature basin paving (MCTBP) method of quantum theory of atoms in molecules (QTAIM) analysis. RSC Advances. 7 (30), 18401-18417 (2017).
  23. Deaven, D. M., Ho, K. M. Molecular geometry optimization with a genetic algorithm. Physical Review Letters. 75, 288-291 (1995).
  24. Hartke, B. Application of evolutionary algorithms to global cluster geometry optimization. Applications of Evolutionary Computation in Chemistry. , Springer. Berlin. (2004).
  25. Dieterich, J. M., Hartke, B. OGOLEM: Global cluster structure optimization for arbitrary mixtures of flexible molecules. A multiscaling, object-oriented approach. Molecular Physics. 108 (3-4), 279-291 (2010).
  26. Herb, J., Nadykto, A. B., Yu, F. Large ternary hydrogen-bonded pre-nucleation clusters in the Earth's atmosphere. Chemical Physics Letters. 518, 7-14 (2011).
  27. Ortega, I. K., et al. From quantum chemical formation free energies to evaporation rates. Atmospheric Chemistry and Physics. 12 (1), 225-235 (2012).
  28. Elm, J., Bilde, M., Mikkelsen, K. V. Influence of Nucleation Precursors on the Reaction Kinetics of Methanol with the OH Radical. Journal of Physical Chemistry A. 117 (30), 6695-6701 (2013).
  29. Loukonen, V., et al. Enhancing effect of dimethylamine in sulfuric acid nucleation in the presence of water - a computational study. Atmospheric Chemistry and Physics. 10 (10), 4961-4974 (2010).
  30. Temelso, B., Phan, T. N., Shields, G. C. Computational study of the hydration of sulfuric acid dimers: implications for acid dissociation and aerosol formation. Journal of Physical Chemistry A. 116 (39), 9745-9758 (2012).
  31. Jiang, S., et al. Study of Cl-(H2O)n (n = 1-4) using basin-hopping method coupled with density functional theory. Journal of Computational Chemistry. 35 (2), 159-165 (2014).
  32. Temelso, B., et al. Effect of mixing ammonia and alkylamines on sulfate aerosol formation. Journal of Physical Chemistry A. 122 (6), 1612-1622 (2018).
  33. Perdew, J. P., Ruzsinszky, A., Tao, J. Prescription for the design and selection of density functional approximations: More constraint satisfaction with fewer fits. Journal of Chemical Physics. 123, 062201(2005).
  34. Riplinger, C., Neese, F. An efficient and near linear scaling pair natural orbital based local coupled cluster method. Journal of Chemical Physics. 138, 034106(2013).
  35. Riplinger, C., Pinski, P., Becker, U., Valeev, E. F., Neese, F. Sparse maps--A systematic infrastructure for reduced-scaling electronic structure methods. II. Linear scaling domain based pair natural orbital coupled cluster theory. Journal of Chemical Physics. 144 (2), 024109(2016).
  36. Kildgaard, J. V., Mikkelsen, K. V., Bilde, M., Elm, J. Hydration of atmospheric molecular clusters: a new method for systematic configurational sampling. Journal of Physical Chemistry A. 122 (22), 5026-5036 (2018).
  37. González, Á Measurement of areas on a sphere Using Fibonacci and latitude-longitude lattices. Mathematical Geosciences. 42, 49-64 (2010).
  38. Karaboga, D., Basturk, B. On the performance of artificial bee colony (ABC) algorithm. Applied Soft Computing. 8 (1), 687-697 (2008).
  39. Zhang, J., Doig, M. Global optimization of rigid molecules using the artificial bee colony algorithm. Physical Chemistry Chemical Physics. 18 (4), 3003-3010 (2016).
  40. Kubecka, J., Besel, V., Kurten, T., Myllys, N., Vehkamaki, H. Configurational sampling of noncovalent (atmospheric) molecular clusters: sulfuric acid and guanidine. Journal of Physical Chemistry A. 123 (28), 6022-6033 (2019).
  41. Grimme, S., Bannwarth, C., Shushkov, P. A Robust and accurate tight-binding quantum chemical method for structures, vibrational frequencies, and noncovalent Interactions of large molecular systems parametrized for all spd-block elements (Z = 1-86). Journal of Chemical Theory and Computation. 13 (5), 1989-2009 (2017).
  42. Buck, U., Pradzynski, C. C., Zeuch, T., Dieterich, J. M., Hartke, B. A size resolved investigation of large water clusters. Physical Chemistry Chemical Physics. 16 (15), 6859(2014).
  43. Forck, R. M., et al. Structural diversity in sodium doped water trimers. Physical Chemistry Chemical Physics. 14 (25), 9054-9057 (2012).
  44. Witt, C., Dieterich, J. M., Hartke, B. Cluster structures influenced by interaction with a surface. Physical Chemistry Chemical Physics. 20 (23), 15661-15670 (2018).
  45. Freitbert, A., Dieterich, J. M., Hartke, B. Exploring self-organization of molecular tether molecules on a gold surface by global structure optimization. The Journal of Computational Chemistry. 40 (22), 1978-1989 (2019).
  46. Stewart, J. J. P. Optimization of parameters for semiempirical methods VI: More modifications to the NDDO approximations and re-optimization of parameters. The Journal of Molecular Modeling. 19 (1), 1-32 (2013).
  47. Stewart, J. J. P. MOPAC2012 Computational Chemistry. , Available from: http://openmopac.net (2012).
  48. Burke, K., Perdew, J. P., Wang, Y. Derivation of a generalized gradient approximation: The PW91 density functional. Electronic Density Functional Theory. , Springer. Boston, MA. 81-111 (1998).
  49. Frisch, M. J., et al. Gaussian 09, Revision A.02. , Gaussian, Inc. Wallingford, CT. (2016).
  50. Ditchfield, R., Hehre, W. J., Pople, J. A. Self-consistent molecular-orbital methods. IX. An extended Gaussian-type basis for molecular-orbital studies of organic molecules. The Journal of Chemical Physics. 54 (2), 724(1971).
  51. Elm, J., Bilde, M., Mikkelsen, K. V. Assessment of density functional theory in predicting structures and free energies of reaction of atmospheric prenucleation clusters. The Journal of Chemical Theory and Computation. 8 (6), 2071-2077 (2012).
  52. Elm, J., Mikkelsen, K. V. Computational approaches for efficiently modelling of small atmospheric clusters. Chemical Physics Letters. 615, 26-29 (2014).
  53. Bayucan, A., et al. PBS Portable Batch System. , MRJ Technology Solutions. Mountain View, CA. (1999).
  54. O'Boyle, N. M., et al. Open Babel: An open chemical toolbox. Journal of Cheminformatics. 3, 33(2011).
  55. Csaszar, A. G. Conformers of gaseous glycine. Journal of the American Chemical Society. 114 (24), 9568-9575 (1992).
  56. Zhang, Q., Anastasio, C. Free and combined amino compounds in atmospheric fine particles (PM2.5) and fog waters from Northern California. Atmospheric Environment. 37 (16), 2247-2258 (2003).
  57. Matsumoto, K., Uematsu, M. Free amino acids in marine aerosols over the western North Pacific Ocean. Atmospheric Environment. 39 (11), 2163-2170 (2005).
  58. Mandalakis, M., Apostolaki, M., Stephanou, E. G. Trace analysis of free and combined amino acids in atmospheric aerosols by gas chromatography-mass spectrometry. Journal of Chromatography A. 1217 (1), 143-150 (2010).
  59. Seinfeld, J. H., Pandis, S. N. Atmospheric Chemistry and Physics, 3rd Ed. , John Wiley & Sons. Hoboken, N.J. (2016).
  60. Myllys, N., Elm, J., Halonen, R., Kurten, T., Vehkamaki, H. Coupled cluster evaluation of atmospheric acid-base clusters with up to 10 molecules. The Journal of Physical Chemistry A. 120 (4), 621-630 (2016).
  61. Elm, J., Bilde, M., Mikkelsen, K. V. Assessment of binding energies of atmospherically relevant clusters. Physical Chemistry Chemical Physics. 15 (39), (2013).
  62. Elstner, M. The SCC-DFTB method and its application to biological systems. Theoretical Chemistry Accounts. 116 (1-3), 316-325 (2006).
  63. Kaliman, I. A., Slipchenko, L. V. LIBEFP: A new parallel implementation of the effective fragment potential method as a portable software library. The Journal of Computational Chemistry. 34 (26), 2284-2292 (2013).
  64. Zhao, Y., Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and trasition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theoretical Chemistry Accounts. 120 (1-3), 215-241 (2008).
  65. Mardirossian, N., Head-Gordon, M. wB97X-V: A 10-parameter, range-separated hybrid, generalized gradient approximation density functional with nonlocal correlation, designed by a survival-of-the-fittest strategy. Physical Chemistry Chemical Physics. 16 (21), 9904-9924 (2014).
  66. Head-Gordon, M., Pople, J. A., Frisch, M. J. MP2 energy evaluation by direct methods. Chemical Physics Letters. 153 (6), 503-506 (1988).
  67. Pople, J. A., Seeger, R., Krishnan, R. Variational configuration interaction methods and comparison with perturbation theory. The International Journal of Quantum Chemistry. 12, 149-163 (1977).
  68. Pople, J. A., Binkley, J. S., Seeger, R. Theoretical models incorporating electron correlation. The International Journal of Quantum Chemistry. 10 (10), 1-19 (1976).
  69. Monkhorst, H. J. Calculation of properties with the coupled-cluster method. The International Journal of Quantum Chemistry. 12 (11), 421-432 (1977).
  70. Klopper, W., Manby, F. R., Ten-No, S., Valeev, E. F. R12 methods in explicitly correlated molecular electronic structure theory. International Reviews in Physical Chemistry. 25, 427-468 (2006).
  71. Hattig, C. Optimization of auxiliary basis sets for RI-MP2 and RI-CC2 calculations: Core-valence and quintuple-z basis sets for H to Ar and QZVPP basis sets for Li to Kr. Physical Chemistry Chemical Physics. 7 (1), 59-66 (2005).
  72. Barone, V. Anharmonic vibrational properties by a fully automated second-order perturbative approach. The Journal of Chemical Physics. 122, 014108(2005).
  73. Barone, V. Vibrational zero-point energies and thermodynamic functions beyond the harmonic approximation. The Journal of Chemical Physics. 120 (7), 3059-3065 (2004).
  74. Temelso, B., et al. Exploring the Rich Potential Energy Surface of (H2O)11 and Its Physical Implications. Journal of Chemical Theory and Computation. 14 (2), 1141-1153 (2018).
  75. Kabrede, H., Hentschke, R. Global minima of water clusters (H2O)N, N≤25, described by three empirical potentials. Journal of Physical Chemistry B. 107 (16), (2003).
  76. Steber, A. L., et al. Capturing the Elusive Water Trimer from the Stepwise Growth of Water on the Surface of a Polycyclic Aromatic Hydrocarbon Acenaphthene. Journal of Physical Chemistry Letters. 8 (23), 5744-5750 (2017).
  77. Perez, C., et al. Corrannulene and its complex with water: A tiny cup of water. Physical Chemistry Chemical Physics. 19 (22), 14214-14223 (2017).

重印与许可

标签

从ab initio热化学计算分子团簇的大气浓度
联系我们
向图书馆推荐

隐私政策

使用条款

政策

版权所有 © 2026 MyJoVE Corporation。保留所有权利。