弱結合分子クラスターの大気濃度は、遺伝的アルゴリズムと半実証的およびab initio量子化学を利用した多段階構成サンプリング方法論を通じて見られる低エネルギー構造の熱化学的特性から計算することができる。
大気エアロゾルの形成と成長の計算研究では、ガス相電子構造と振動周波数計算から得ることができる正確なギブス自由エネルギー表面が必要です。これらの量は、ジオメトリが潜在的なエネルギーサーフェスの最小値に対応する大気クラスターに対して有効です。最小エネルギー構造のギブス自由エネルギーは温度および圧力のようないろいろな条件の下でクラスターの大気の集中を予測するために使用することができる。遺伝的アルゴリズムベースの構成サンプリングに基づいて構築された計算上の安価な手順を提示し、その後で一連のますます正確なスクリーニング計算を行います。この手順は、半経験的モデルを使用して大規模な構成セットの幾何学を生成および進化させることから始まり、結果として生じる一意の構造を一連の高レベルの ab initio理論レベルで洗練します。最後に、熱力学的補正は、生成された最小エネルギー構造のセットに対して計算され、ギブズの形成、平衡定数、大気濃度の自由エネルギーを計算するために使用されます。この手順を、周囲条件下での水和グリシンクラスターの研究に応用して紹介する。
気候変動の大気研究における最も不確実なパラメータは、雲粒子が入ってくる日射量をどの程度反射させるかです。気体に懸濁した粒子状物質であるエアロゾルは、入射放射線を散乱する雲結縮核(CCN)と呼ばれる雲粒子を形成し、大気1の吸収およびその後の加熱を防止する。この正味冷却効果を詳細に理解するには、エアロゾルのCCNへの成長を理解する必要があり、その結果、小分子クラスターがエアロゾル粒子に成長するのを理解する必要があります。最近の研究では、エアロゾル形成は直径3nm以下の分子クラスターによって開始される2;しかし、このサイズレジームは、実験技術33、44を使用してアクセスすることは困難です。したがって、この実験上の限界を克服するために、計算モデリング手法が望まれている。
以下に説明するモデリング手法を使用して、水和されたクラスターの成長を分析できます。我々は、生物前環境における小さな構成物質からの大きな生体分子の形成における水の役割に関心を持っているので、グリシンによるアプローチを示す。これらの研究課題に対処するために必要な課題とツールは、大気エアロゾルと前,核クラスタ,55、6、7、8、9、10、11、12、13、14、156,7,の研究に関与するものと非常によく似ています。8,910,1112,13,14,15ここでは、単離されたグリシン分子から始まる水和グリシンクラスターを調べ、続いて最大5つの水分子を段階的に添加した。最終目標は、大気圏のGly(H2O)n=0-5クラスターのn=0-5平衡濃度を、海面で室温、相対湿度(RH)100%を計算することです。2
これらのサブナノメートル分子クラスターの少数は、他の蒸気分子を追加するか、既存のクラスター上で凝固することによって、メタスタブルクリティカルクラスター(直径1〜3nm)に成長します。これらの重要なクラスターは、はるかに大きな(最大50〜100 nm)雲凝縮核(CCN)の形成につながる良好な成長プロファイルを有し、雲の降水効率に直接影響を及ぼし、入射光を反射する能力も直接影響を及ぼします。したがって、分子クラスターの熱力学とその平衡分布を理解することは、エアロゾルが地球規模の気候に与える影響をより正確に予測する必要があります。
エアロゾル形成の記述モデルには、分子クラスター形成の正確な熱力学が必要です。分子クラスター形成の正確な熱力学の計算では、最も安定した構成の同定が必要であり、これは、クラスターの潜在的エネルギー表面(PES)16上のグローバルおよび局所的なミニマを見つけることを含む。このプロセスは、構成サンプリングと呼ばれ、分子動力学(MD)17、18、19、20、モンテ,20カルロ18,19(MC)21、22、および遺伝的21,アルゴリズム17,(GA)23、24、25に基づくもの23を含む様々な技術,24を通じて達成することができる。2225
高い理論レベルで大気水和物の構造と熱力学を得るために、長年にわたって異なるプロトコルが開発されてきました。これらのプロトコルは、(i)構成サンプリング法、(ii)構成サンプリングで使用される低レベル法の性質、および(iii)後のステップで結果を洗練するために使用される高レベルのメソッドの階層の選択が異なった。
構成サンプリング法には、化学的直感26、ランダムサンプリング27、28、28分子動力学(MD)29、30、洗面体ホ29,30ッピング(BH)31、および遺伝的アルゴリズム(GA)24、25、32が含まれていた32。24,25,2731これらのサンプリング方法で使用される最も一般的な低レベルの方法は、力場またはPM6、PM7およびSCC-DFTBのような半経験的モデルである。これらはしばしば、ヤコブのはしご33の高いラングからますます大きなベーシスセットとより信頼性の高い機能を持つDFT計算が続きます。場合によっては、MP2、CCSD(T)、コスト効率の高い DLPNO-CCSD(T)34、,35など、より高いレベルの波動関数方式が続きます。
Kildgaardら36は、フィボナッチ球体37の小さな水和または非水和クラスターの周りに水分子を添加し、より大きなクラスターの候補を生成する体系的な方法を開発した。非物理的な候補と冗長な候補は、近接接触しきい値と異なる適合者間のルート平均平方距離に基づいて削除されます。PM6セミ経験的手法とDFTと波関数法の階層を用いた後の最適化は、高い理論レベルで低エネルギー適合者のセットを得るために使用される。
人工的なハチコロニー(ABC)アルゴリズム38は、最近、ABCluster39と呼ばれるプログラムで分子クラスターを研究するためにZhangらによって実施された新しい構成サンプリングアプローチである。Kubeckaら40は、構成サンプリングのためにABClusterを使用し、続いて、緊密結合GFN-xTB半実証法41を用いた低レベルの再最適化を行った。さらに、DFT法を用いて構造とエネルギーを改良し、その後にDLPNO-CCSD(T)を用いた最終的なエネルギーを使用しました。
方法に関係なく、構成サンプリングは、ペスのポイントのランダムまたは非ランダムに生成された分布から始まります。各点は、問題の分子クラスターの特定の幾何学に対応し、サンプリング法によって生成されます。次に、PESの「下り坂」方向に従うことによって、各ポイントに最も近い局所最小値が見つかります。このように見つかったミニマのセットは、少なくともしばらくの間、安定している分子クラスターの幾何学に対応する。ここで、PESの形状と、表面上の各点におけるエネルギーの評価は、より正確な物理的記述がより計算的に高価なエネルギー計算をもたらすシステムの物理的記述に敏感であろう。我々は、特に、サンプリングポイントの初期セットを生成するために25、グローバルな最適化と構成サンプリングの問題42、43、44、4543,44,の42様々なに適用されているOGOLEM 25プログラムで実装されたGAメソッドを使用します。45PESは、MOPAC2016プログラム47に実装されたPM7モデル46で説明する。この組み合わせは、MDおよびMC法に比べてより多様なポイントを生成し、PESのより詳細な記述よりもローカルミニマを迅速に検出するため採用されています。
GA 最適化されたローカル ミニマのセットは、一連のスクリーニングステップの開始ジオメトリとして取られ、低位置の最小エネルギーのセットにつながります。プロトコルのこの部分は、密度機能理論(DFT)を使用して、小さな基本セットを使用して、ユニークなGA最適化構造のセットを最適化することから始まります。この最適化のセットは、一般的に、GA最適化半経験的構造と比較して、より詳細にモデル化される一意のローカル最小構造の小さなセットを与えます。次に、より大きな基準セットを使用して、この小さな構造体のセットに対してDFT最適化の別のラウンドが実行されます。繰り返しますが、このステップは一般的に、小さな基礎DFTステップと比較してより詳細にモデル化されるユニークな構造の小さなセットを与えます。その後、固有の構造の最終的なセットは、より緊密な収束に最適化され、調和振動周波数が計算されます。このステップの後、我々は大気中のクラスターの平衡濃度を計算するために必要なすべてを持っています。全体的なアプローチは図1に図示されています。ガウシアン49DFT実装におけるPW9148一般化勾配近似(GGA)交換相関関数と、Pople50ベーシスセット(小ベーシスステップでは6-31+G*、大ベーシスステップでは6-311++G**)の2つのバリエーションを使用します。交換相関機能とベーシスセットのこの特定の組み合わせは、大気クラスター51、52,52のための形成の正確なギブス自由エネルギーを計算する上での以前の成功のために選ばれました。
このプロトコルは、ポータブル バッチ システム53 (PBS)、MOPAC2016 (http://openmopac.net/MOPAC2016.html) 47、OGOLEM (https://www.ogolem.org)25、ガウス 09 (https://gaussian.com)49、および特定のインストール49手順に従ってインストールされた OpenBabel54 (http://openbabel.org/wiki/Main_Page) ソフトウェアを使用して、高性能コンピューティング (HPC) クラスターにアクセスできるものとします。このプロトコルの各ステップでは、社内シェルと Python 2.7 スクリプトのセットも使用しており、ユーザーの$PATH環境変数に含まれるディレクトリに保存する必要があります。上記のすべてのプログラムを実行するために必要なすべての環境モジュールと実行権限も、ユーザーのセッションにロードする必要があります。GA コード (OGOLEM) と準経験コード (MOPAC) によるディスクとメモリの使用量は、最新のコンピュータ リソース標準では非常に小さいものです。OGOLEM/MOPAC の全体的なメモリとディスク使用量は、使用するスレッドの数によって異なり、ほとんどの HPC システムの機能に比べてリソースの使用量は小さくなります。QM法のリソースニーズは、クラスターのサイズと使用される理論のレベルによって異なります。このプロトコルを使用する利点は、理論のレベルを変えて、低エネルギー構造の最終セットを計算できることです。
ユーザーのローカル コンピュータは、わかりやすくするため、”ローカル コンピュータ” と呼ばれ、アクセスできる HPC クラスタは “リモート クラスタ” と呼ばれます。
このプロトコルによって生成されるデータの精度は、主に、(i)ステップ2でサンプリングされた構成の多様性、(ii)システムの電子構造の精度、(iii)および熱力学的補正の精度の3つの要素に依存する。これらの各要因は、含まれているスクリプトを編集してメソッドを変更することで対処できます。最初の要因は、ランダムに生成された構造のより大きな初期プール、GAのより多くの反復、およびGAに関与する基準のより緩い定義を使用して簡単に克服されます。また、異なる物理的記述の影響を探るために、自己一貫した電荷密度機能密着性の密結合(SCC-DFTB)6262モデルおよび有効断片電位(EFP)6363モデルのような異なる半経験的方法を用いてもよい。ここでの主な制限は、共有結合を形成または切断する方法の不能であり、単量体が凍結していることを意味する。GA手順は、半経験的記述に従って、これらの凍結モノマーの最も安定した相対位置のみを見つける。
システムの電子構造の精度は、その計算コストを持つ様々な方法で改善することができます。M06-2X 64やwB97X-V65などのより優れた密度機能、65またはMøller-Plesset69 66、67、68(MPn)66摂動理論などの量子機械(QM)法を選択して、67,68システムの物理的記述を改善することができます。機能の階層では、一般的に、PW91のような一般化された勾配近似(GGA)機能から、M06-2XのようなwB97X-DおよびメタGGAハイブリッド機能のような範囲分離されたハイブリッド機能に至ると、性能が向上する。
DFT法の欠点は、正確な値に向けた体系的な収束が不可能であるということです。しかし、DFT法は計算上安価であり、多種多様なアプリケーションに対して多種多様な機能があります。
増加する枢機卿数の相関一貫基準セットと組み合わせてMP2やCCSD(T)のような波関数法を使用して計算されたエネルギー([aug-]cc-pV[D,T,Q,…]Z) は、完全な基準セットの制限に向かって体系的に収束しますが、システムサイズが大きくなるにつれて、各計算の計算コストは非常に大きくなります。電子構造のさらなる改良は、明示的に相関するベーシスセット70を用いて、完全な基準セット(CBS)71制限に外挿することによって達成71することができる。我々の最近の研究は、密度が明示的に相関する第2次Møller-Plesset(DF-MP2-F12)摂動アプローチがMP2/CBS計算32のエネルギーに近づくエネルギーをもたらすことを示唆している。現在のプロトコルを変更して異なる電子構造メソッドを使用するには、(i)ソフトウェアの構文に従ってテンプレート入力ファイルを作成し、(ii)run-pw91-sb.csh、run-pw91-lb.csh、およびrun-pw91-lb-ultrafine.cshスクリプトを使用して正しい入力ファイル構文を生成する2つの手順を実行します。 run-pw91-lb.csh
最後に、熱力学による補正の精度は、電子構造法と、グローバル最小の周りのPESの記述に依存します。PESの正確な記述では、非常にコストのかかる作業である、四次力場72,73(QFF),73のような核自由度の変位に関するPESの3次および高次微分の計算が必要である。現在のプロトコルは振動周波数に対して高調波発振器近似を使用するため、PESの二次導関数までしか計算する必要が出ない。このアプローチは、真のPESと高調波PESの大きな違いのために、非常にフロッピー分子や対称的なダブルウェル電位など、高い不均調性を持つシステムでは問題になります。さらに、計算負荷の高い電子構造法から高品質PESを有するコストは、振動周波数計算のコストの問題を複合化するのみである。これを克服する1つのアプローチは、高品質の電子構造計算から電子エネルギーを使用し、低品質のPESで計算された振動周波数と共に、コストと精度のバランスを取ることです。現在のプロトコルは、前の段落で説明したように、異なるPESの記述を使用するように変更することができます。ただし、スクリプトやテンプレートの振動周波数キーワードを編集して、アンコン波振動周波数を計算することもできます。
構成サンプリング プロトコルの 2 つの重要な問題は、潜在的なエネルギー サーフェスをサンプリングするための最初の方法と、各クラスターを識別するために使用される基準です。前作では様々な方法を幅広く活用してきました。最初の問題として、潜在的なエネルギー表面をサンプリングするための最初の方法は、これらの要因に基づいて半経験的な方法でGAを使用するという選択をしました。化学直感26、ランダムサンプリング、分子動力学(MD)29,30を用いた構成サンプリングは、水29,30クラスター18の研究で観察したように、10個のモノマーを超えるクラスターに対して定期的に推定グローバルミニマを見つけることができません。我々は、(H2O)111174の複雑なPESを研究するために2盆地ホッピング(BH)を使用することに成功したが、それはBHアルゴリズムが見つけ出さなかったいくつかの潜在的な低エネルギー異性体の手動包含を必要とした。水クラスターの世界的な最小値を見つける場合のBHとGAのパフォーマンスの比較, (H2O)n=10-20 GAは一貫してBH75よりも速いグローバル最小値を発見したことを示しました.OGOLEM および CLUSTER に実装されている GA は、任意の分子クラスターに適用でき、古典的な力場、半経験的、密度機能、および ab initio機能を持つ膨大な数のパッケージとのインターフェースが可能であるため、非常に汎用性があります。PM7の選択は速度および合理的な正確さによって動かされる。他のセミ経験的手法は、計算コストが大幅に高くなります。
2つ目の課題としては、電子エネルギー、双極子モーメント、オーバーラップRMSD、回転定数に至るまで、異なる基準を用いて、独自の構造を特定しました。双極子モーメントの両方が分子の向きに依存し、ダイポールモーメントの総部分が、構造が同じまたはユニークであると判断する閾値を設定することが困難であったような形状の違いに非常に敏感であったため、双極子モーメントを使用することは困難であることが判明した。電子エネルギーと回転定数の組み合わせが最も有用であることが判明しました。
2つの構造を一意と見なす現在の基準は、エネルギー差閾値0.10kcal mol-1と回転定数差1%に基づいています。したがって、エネルギーが0.10kcalモル-1(〜0.00015-1 a.u.)以上異なる場合、2つの構造が異なると考えられます。そして、その3つの回転定数(A、B、C)のいずれかが1%以上異なります。長年にわたる実質的な内部ベンチマークは、これらのしきい値が合理的な選択であることを発見しました。我々の構成サンプリング手法及びスクリーニング方法論は、水76,77,77と、アンモニア及びアミン類を含む強結合した三硫酸水和物と同様に、水と結合した多芳香族炭化水素のような非常に弱く結合したクラスターに適用された。考慮するプロネーション状態が異なるクラスターの場合、各々が異なる原始状態のモノマーから始まるさまざまな GA 計算を実行するのが最善の方法です。これにより、異なる原始状態の構造を慎重に考慮することが保証されます。ただし、低レベル DFT 計算では、多くの場合、ジオメトリ最適化の過程でプロトネーション状態が変化するため、開始ジオメトリに関係なく最も安定したプロトネーション状態が得られます。
GA コードが一般的な非パラメータ化メソッドとインターフェイスされている限り、GA の GA 構成サンプリング方法は、GA の実行中にモノマーが異なる構成を採用できる限り、うまく機能するはずです。例えば、GAとPM7をインターフェースするとモノマーの構造が変わることになりますが、原始状態が変化したときに結合が壊れた場合、構造は受け入れられない候補として廃棄される可能性があります。
高調波近似の欠点、特に低振動周波数に起因する欠点を修正するさまざまな方法を検討してきました。準高調波近似を現在の方法論に組み込むことは難しくない。しかし、準高調波法については、特に適用されるカットオフ周波数に関しては疑問が残っています。また、従来の知恵がRRHO近似値よりも改善すべきであることを示唆しているにもかかわらず、準RRHO近似の信頼性を調べる厳格なベンチマーク作業はありません。
このように提示されるプロトコルは、非共有結合ガス相分子クラスターの任意のシステムに一般化され得る。また、スクリプトやテンプレートを編集することにより、セミ経験的方法、電子構造法およびソフトウェア、および振動解析方法およびソフトウェアを使用するように一般化してもよい。これは、ユーザーが Linux コマンドライン インターフェイス、Python スクリプト、およびハイパフォーマンス コンピューティングに慣れかけられていることを前提としています。Linux オペレーティングシステムの不慣れな構文と見た目とスクリプトの経験の欠如は、このプロトコルの最大の落とし穴であり、新入生が最も苦労している場所です。このプロトコルは、主に硫酸とアンモニアがエアロゾル形成に及ぼす影響に焦点を当てて、私たちのグループで長年にわたって様々な実装で成功してきました。このプロトコルのさらなる改善には、より多くの電子構造ソフトウェアへのより堅牢なインターフェース、遺伝的アルゴリズムの代替実装、そしておそらく電子および振動エネルギーのより速い計算のための新しい方法の使用が含まれる。このプロトコルの現在のアプリケーションは、現在の大気中でのエアロゾル形成の初期段階とプレバイオティクス環境におけるより大きな生体分子の形成におけるアミノ酸の重要性を探っています。
The authors have nothing to disclose.
このプロジェクトは、国立科学財団(GCS)、アーノルドとメイベルベックマン財団ベックマン奨学生賞(AGG)、バリーM.ゴールドウォーター奨学金(AGG)からの助成金CHE-1229354、CHE-1662030、CHE-1721511、CHE-1903871によって支援されました。水銀コンソーシアム(http://www.mercuryconsortium.org)の高性能コンピューティング資源を使用した。
Avogadro | https://avogadro.cc | Open-source molecular visualization program | |
Gaussian [09/16] Software | http://www.gaussian.com/ | Commercial ab initio electronic structure program | |
MOPAC 2016 | http://openmopac.net/MOPAC2016.html | Open-source semi-empirical program | |
OGOLEM Software | https://www.ogolem.org | Genetic algorithm-based global optimization program | |
OpenBabel | http://openbabel.org/wiki/Main_Page | Open-source cheminformatics library | |
calcRotConsts.py | Shields Group, Department of Chemistry, Furman University | Python script to compute rotational constants | |
calcSymmetry.csh | Shields Group, Department of Chemistry, Furman University | Shell script to calculate symmetry number of a molecule given Cartesian coordinates | |
combine-GA.csh | Shields Group, Department of Chemistry, Furman University | Shell script to combine energy and rotational constants from different GA directories | |
combine-QM.csh | Shields Group, Department of Chemistry, Furman University | Shell script to combine energy and rotational constants from different QM directories | |
gaussianE.csh | Shields Group, Department of Chemistry, Furman University | Shell script to extract Gaussian 09 energies | |
gaussianFreqs.csh | Shields Group, Department of Chemistry, Furman University | Shell script to extract Gaussian 09 vibrational frequencies | |
getrotconsts | Shields Group, Department of Chemistry, Furman University | Executable to calculate rotational constants given a molecule's Cartesian coordinates | |
getRotConsts-dft-lb.csh | Shields Group, Department of Chemistry, Furman University | Shell script to compute rotational constants for a batch of large basis DFT optimized structures | |
getRotConsts-dft-lb-ultrafine.csh | Shields Group, Department of Chemistry, Furman University | Shell script to compute rotational constants for a batch of ultrafine DFT optimized structures | |
getRotConsts-dft-sb.csh | Shields Group, Department of Chemistry, Furman University | Shell script to compute rotational constants for a batch of small basis DFT optimized structures | |
getRotConsts-GA.csh | Shields Group, Department of Chemistry, Furman University | Shell script to compute rotational constants for a batch of genetic algorithm optimized structures | |
global-minimum-coords.xyz | Shields Group, Department of Chemistry, Furman University | Cartesian coordinates of global minimum structures of gly-(h2o)n, where n=0-5 | |
make-thermo-gaussian.csh | Shields Group, Department of Chemistry, Furman University | Shell script to extract data from Gaussian output files and make input files for the thermo.pl script | |
ogolem-input-file.ogo | Shields Group, Department of Chemistry, Furman University | Ogolem sample input file | |
ogolem-submit-script.pbs | Shields Group, Department of Chemistry, Furman University | PBS batch submission file for Ogolem calculations | |
README.docx | Shields Group, Department of Chemistry, Furman University | Clarifications to help readers use the scripts effectively | |
runogolem.csh | Shields Group, Department of Chemistry, Furman University | Shell script to run OGOLEM | |
run-pw91-lb.csh | Shields Group, Department of Chemistry, Furman University | Shell script to run a batch of large basis DFT optimization calculations | |
run-pw91-lb-ultrafine.csh | Shields Group, Department of Chemistry, Furman University | Shell script to run a batch of ultrafine DFT optimization calculations | |
run-pw91-sb.csh | Shields Group, Department of Chemistry, Furman University | Shell script to run a batch of small basis DFT optimization calculations | |
run-thermo-pw91.csh | Shields Group, Department of Chemistry, Furman University | Shell script to compute the thermodynamic corrections for a batch of DFT optimized structures | |
similarityAnalysis.py | Shields Group, Department of Chemistry, Furman University | Python script to determine unique structures based on rotational constants and energies | |
symmetry | Shields Group, Department of Chemistry, Furman University | Executable to calculate molecular symmetry given Cartesian coordinates | |
symmetry.c | (C) 1996, 2003 S. Patchkovskii, Serguei.Patchkovskii@sympatico.ca | C code to determine the molecular symmstry of a molecule given Cartesian coordinates | |
template-marcy.pbs | Shields Group, Department of Chemistry, Furman University | Template for a PBS submit script which uses OGOLEM | |
template-pw91.com | Shields Group, Department of Chemistry, Furman University | Template Gaussian 09 input | |
template-pw91-HL.com | Shields Group, Department of Chemistry, Furman University | Template Gaussian 09 input for ultrafine DFT optimization | |
thermo.pl | https://www.nist.gov/mml/csd/chemical-informatics-research-group/products-and-services/program-computing-ideal-gas | Perl open-source script to compute ideal gas thermodynamic corrections | |
gly-h2o-n.xlsx | Shields Group, Department of Chemistry, Furman University | Excel spreadsheet for the complete protocol | |
table-1.xlsx | Shields Group, Department of Chemistry, Furman University | Excel spreadsheet | |
table-2.xlsx | Shields Group, Department of Chemistry, Furman University | Excel spreadsheet | |
table-3.xlsx | Shields Group, Department of Chemistry, Furman University | Excel spreadsheet | |
table-4.xlsx | Shields Group, Department of Chemistry, Furman University | Excel spreadsheet | |
water.xyz | Shields Group, Department of Chemistry, Furman University | Cartesian coordinates of water | |
glycine.xyz | Shields Group, Department of Chemistry, Furman University | Cartesian coordinates of glycine |