<100>引張荷重下のタンタル単結晶の弾性特性の温度および圧力依存性：分子動力学研究

背景

一般に、タンタル（Ta）は周囲条件でBCC構造に属します。現在、多くの文献が、Ta単結晶が高圧下で優れた相安定性を示すことを証明しています[1,2,3]。さらに、Taの周囲圧力での溶融温度は3269 Kと非常に高く、他のほとんどの金属よりも高くなっています[4]。 Taはその優れた特性により、マイクロ/ナノエレクトロニクスの拡散バリア、耐摩耗性コーティング、高温超合金など、多くの技術アプリケーションに理想的な材料です。

最近、実験[2、5、6、7]と理論[8、9、10、11、12、13、14]の両方の分野で多大な努力が、 Ta。 Dewaele etal。 [5]は、93 GPaまでのダイヤモンドアンビルセル（DAC）のTaの降伏強度に対する圧力の影響を研究し、DAC実験では、BCC構造が135GPaまで安定していることも示されました[2]。さらに、Shigeaki [8]は、DFTを使用して最大100GPaおよび3000KのTaの状態方程式（EOS）をシミュレートしました。ウーら[9]は、350GPaまでの高圧でのTaの弾性および熱力学的特性を調査しました。一方、Škoroetal。 [6、7]は、それぞれ2250Kと2500Kまでの非常に高い温度でのTaの降伏強度とヤング率を測定しました。 Gu etal。 [10]は、第一原理法を使用して、500GPaまでの立方晶Taの高圧構造と弾性特性に関する研究を実施しました。圧力の関数としての弾性定数、およびTaのバルク、ヤング率、およびせん断弾性率はすべて、圧力の増加とともに増加することがわかります。

DAC実験とDFT計算に加えて、MDシミュレーションの分野での高温高圧での多くの研究もあります[15、16、17、18]。 Liu etal。 [15]は、拡張されたFinnis-Sinclair（EFS）ポテンシャルを使用し、最大400GPaの圧力でのTaの熱EOSと溶融特性を調査しました。その上、Tramontina等。 [16、17]は、高圧での可塑性メカニズムに対する結晶配向の影響を研究しました。彼らはまた、衝撃強度と衝撃立ち上がり時間がそれらの微細構造に及ぼす影響についても議論しました。さらに、Ruestes etal。 [18] 3つの異なる原子間ポテンシャルを使用してBCCTaの押し込みシミュレーションを実行し、塑性変形ゾーンの作成と拡張の原因となる欠陥メカニズムを示しました。

上記の多くの調査にもかかわらず、MDシミュレーションを使用した引張荷重下でのTaの動的応答の体系的な原子シミュレーション研究はありませんでした。本研究の主な目的は、サイズ、ひずみ速度、温度、および圧力の影響を考慮して、<100>引張荷重下でのTa単結晶の弾性特性を調査することです。さらに、<100>引張荷重によって相転移が誘発されるかどうかを理解することも、この作業のもう1つの目的です。

メソッド/実験

物理モデリング

図1に示すように、この論文で調査したTaキュービックは、<001>、<010>、および<100>の方向に沿ってBCCユニットセルを繰り返すことによって生成され、格子定数は a です。 = b = c =3.301Å、それぞれ。 12 a を含む、エッジの長さが異なる4つの立方体モデル （3.96 nm）、18 a （5.94 nm）、24 a （7.92 nm）、および30 a （9.90 nm）、構築されます。対応する原子数は、それぞれ3456、11,664、27,648、および54,000です。図1は、エッジ長が3.96 nmのTa立方体のスケッチマップを示しています。これは、現在のシミュレーションの元の構造です。

エッジ長が3.96nmのTa立方体のBCC構造とスケッチマップ。スケッチマップは、現在のシミュレーションの元の構造です

MDシミュレーションの詳細

この論文で報告されているMDシミュレーションは、大規模な原子/分子の大規模並列シミュレータ（LAMMPS）を使用して実行されます[19]。 MDシミュレーションの過程で、Ta原子に加えられる外力は、それらの原子間の原子間ポテンシャル関数に従って計算されます。この研究では、2つの異なる原子間ポテンシャルが考慮されています。ZhouらによるEAMポテンシャル。 [20]およびRaveloetal。それぞれ[21]。両方のポテンシャルは、表1にリストされているTaのいくつかの弾性特性を再現するために適用されました。簡単にするために、これらはRavelo-EAMおよびZhou-EAMポテンシャルと呼ばれます。

立方対称の構造の弾性挙動は、弾性定数 C によって完全に記述されます。 ₁₁ 、 C ₁₂ 、および C ₄₄ 。 <100>、<110>、および<111>単結晶の配向依存弾性率は、いくつかの方程式[18]を使用して計算されます。私たちの研究では、MDシミュレーションを実行して、<100>引張荷重下での弾性率への影響を調べています。したがって、<100>方向の弾性率に焦点を当て、 E の弾性率を採用します。 ₁₀₀ 考慮に入れます。したがって、弾性定数 C を考慮します。 ₁₁ および C ₁₂ また、次の方程式[22]：

一般に、MDから弾性定数を計算するには、応力変動法、ひずみ変動法、直接法の3つの方法があります。現在の研究では、Gaoらの方法と同様の直接法。 [23]は、弾性定数 C の計算に使用されます。 ₁₁ および C ₁₂ 図2に示すように、2種類のEAMポテンシャルの場合。

C の弾性定数の圧力依存性 ₁₁ および C ₁₂ 研究された可能性に基づいています。 a 周EAM。 b Ravelo-EAM。本研究では、直接法を利用して弾性定数 C を計算しました。 ₁₁ および C ₁₂ 2種類のEAMポテンシャル。 2つの電位のC11 / C12-圧力曲線を a に示します。 および b 、それぞれ

図2aから、両方の C の得られた曲線を見つけることができます。 ₁₁ および C ₁₂ Zhou-EAMの可能性に基づくと、140GPaまでの圧力の増加に伴ってスムーズに進むことはできません。 〜40 GPaの圧力でクロスオーバーポイントがあります。つまり、 C ₁₁ = C ₁₂ 、 E の弾性率 ₁₀₀ 式に従ってゼロになります。 （1）。さらに、 E ₁₀₀ 40 GPaを超える圧力では負の値が示されますが、これは疑わしく、理論的および実験的な結果と矛盾します[24、25]。したがって、Zhou-EAMポテンシャルは、ここで検討した範囲では十分に機能しません。次に、 C のMDシミュレーション結果から、Ravelo-EAMポテンシャルの実現可能性を調べてみましょう。 ₁₁ および C ₁₂ 図2bに示されているRavelo-EAMポテンシャルに基づいています。数値結果は、圧力が高いほど、両方の C の値が大きくなることを示しています。 ₁₁ および C ₁₂ 、これは、DFT計算による圧力によって変化する弾性定数の傾向とよく一致します[9、25]。一方、Ravelo-EAMポテンシャルを使用して計算された結果は、Ruestes etalによって報告された値と同じように注目に値します。 [18]。 Ravelo-EAMポテンシャルは、高圧下で良好に動作すると同時に、動的変形下でのTaの弾性および機械的特性を説明することもできます[26]。したがって、次のセクションでは、Ravelo-EAMの可能性に基づいてシミュレーションを実行します。

幾何学的構築の後、一連の関連するMDシミュレーションを実行します。 MDシミュレーションでは、周期境界条件（PBC）が3次モデルの3方向すべてで使用されます。 Taの弾性特性の温度依存性を調べるために、時間ステップを1 fsに設定し、システム温度を1、300、600、900、1200、および1500Kに設定します。まず、システムが局所的なポテンシャルの最小値にあるように、正準集団（NVT）MDシミュレーションを使用してモデルを約50psの緩和プロセスで緩和します。次に、等温定圧（NPT）MDシミュレーションを使用して、0〜140 GPaの範囲の指定された静水圧を確保し、Taの弾性特性に対する圧力の影響を調べました[27]。最後に、ひずみ速度が5×10 ⁸ の範囲の引張荷重 s ^{− 1} 〜7.5×10 ⁹ s ^{− 1} [28、29]はTaキュービックのx方向に適用されます。一方、NPTシミュレーションは、2番目のステップで適用された同じ圧力でy方向とz方向に実行されます。したがって、ここで計算される弾性率は、<100>方向の場合です。すべてのMDシミュレーションで、モデルは<100>引張荷重によってx方向に15％の伸びまで引き伸ばされます。

結果と考察

ストレッチプロセス

ストレッチプロセス中に、シミュレートされた構成は、科学ソフトウェアパッケージのOpen Visualization Tool（OVITO）[30]を使用して視覚化されます。ゼロ圧力での<100>一軸引張ひずみ下でのTaの応力-ひずみ曲線と、異なるひずみでの対応する原子構成を図3に示します。

ゼロ圧力でのTaの応力-ひずみ曲線と、異なるひずみでの対応する構成。注：青、緑、白のボールは、それぞれBCC、FCC、その他の格子構造に対応しています

図3に示すように、応力-ひずみ曲線から、構成（IV）の近くで表面が破壊されていると結論付けることができます。一軸引張ひずみの開始時に、応力は通常、ひずみに比例して変化します。図3Iは、原子構成がBCC構造を維持していることを示しています。図3II、IIに示すように、ひずみが増加すると、BCCからFCC構造への相転移は、それぞれ〜7.4％のひずみで始まり、〜9.8％のひずみで完了します。そして、これらのFCC構造は、最初の表面破壊の前に最大値を維持します。ひずみが〜13.1％の場合、y方向とz方向のエッジ長が急激に減少し、表面破壊が発生します。一方、図3IVに示すように、クラスターが非常に短時間で発生することは注目に値します。連続的な一軸変形下では、ひずみが約13.3％になるまで原子配置はストライプ状に維持されます。これを図3Vに示します。

この論文では、このセクションで説明するように、引張特性のモデルサイズ、ひずみ速度、温度、および圧力への依存性に焦点を当てます。理論的には、モデルは弾性変形段階で線形に引き伸ばされ、弾性係数は応力-ひずみ曲線の線形部分の勾配として定義されます。すべてのモデルが同様の伸縮プロセスを持ち、応力-ひずみ曲線が同様の​​変動傾向を持っていることがわかります。したがって、同じアプローチを使用して、さまざまなモデルサイズとひずみ速度のTaの弾性率を取得します。

サイズとひずみ速度への依存

表2に、温度1Kおよびひずみ速度5×10 ⁸ でのさまざまなモデルサイズの弾性率と降伏応力を示します。 s ^{− 1} 。モデルサイズは、Taの弾性率と降伏応力に影響を与えないと簡単に結論付けることができます。これは、弾性率が原子間の相互作用を表すことであり、弾性率はモデルのサイズによって変化しないことを説明するのは非常に簡単です。表2から、弾性率は〜139 GPaであり、140GPaのシミュレーション結果と著しく同じであることがわかります[18]。

既存のレフリーによると。 [28、29]、ほとんどのひずみ速度は10 ⁸ の範囲です。 s ^{− 1} 〜10 ¹⁰ s ^{− 1} 。この論文では、5.0×10 ⁸ を含む4つのひずみ速度を選択して引張シミュレーションを実行します。 s ^{− 1} 、7.5×10 ⁸ s ^{− 1} 、5.0×10 ⁹ s ^{− 1} 、および7.5×10 ⁹ s ^{− 1} 。表3に、300Kの温度とさまざまなひずみ速度での弾性率と降伏応力を示します。ひずみ速度は弾性率と降伏応力に明らかな影響を与えないと簡単に結論付けることができます。

一方、さまざまな温度と圧力下での弾性率と降伏応力に対するモデルサイズとひずみ速度の影響もシミュレートします。これらのシミュレーションは同じ結論に達します。したがって、3.96nmの同じモデルサイズと5×10 ⁸ の同じひずみ速度を使用します。 s ^{− 1} 次のシミュレーション用。

温度への依存

図4は、1500 Kまでのさまざまな温度での応力-ひずみ曲線を示しています。これらの曲線の勾配は、<100>方向、つまり E の弾性率を示していることがわかります。 ₁₀₀ 、弾性引張期間中、降伏応力は温度の上昇とともに徐々に減少します。熱力学の理論[31]によると、システムのすべての原子の総運動エネルギーは、一般に次の方程式を満たします。

ここで E _k システムの総運動エネルギーです。 N は原子の総数です。 K _B ボルツマン定数です。 T は熱力学的温度です。したがって、システムはより高い温度でより大きな総運動エネルギーを含み、原子はより速く移動すると結論付けることができます。熱力学的観点から、原子はより活発になり、原子の動きはより激しくなります。これは、平衡位置での振幅が大きくなることを意味します。伸縮過程では、原子間の引力が比較的小さくなり、原子が平衡位置から逃げやすくなるため、同じひずみでx方向の応力が小さくなります。したがって、高温での弾性率は低温での弾性率よりも小さくなります。さらに、これらの曲線の傾向は、Taの初期の理論的および実験的発見とよく一致しています[6、7、32]。

5×10 ⁸ のひずみ速度での3.96nmTaモデルの応力-ひずみ曲線 s ^{− 1} 図4から、これらの曲線の傾きは、<100>方向、つまり E の弾性率を示していることがわかります。 ₁₀₀ 、弾性引張期間中、降伏応力は温度の上昇とともに徐々に減少します

観察を容易にするために、表4に、さまざまな温度でのTaの弾性率と降伏応力を示します。弾性率は136.49から76.67GPaに〜42.3％減少し、降伏応力は温度が1から1500 Kに上昇するにつれて〜8から〜4 GPaに〜51％減少します。

表4によると、弾性率（ E ）の温度依存性のパラメーター化をさらに導き出すことができます。 ₁₀₀ ）結果は図5aに実線で示されています。パラメータ化は次のようになります。

ここで E ₁₀₀ （GPa）と T にあります （ K ）; a =138.07±0.92111、および b =−0.04094±0.00101。この式は、 E ₁₀₀ 温度の上昇とともに直線的に減少するため、0〜1500Kの温度まで使用することをお勧めします。 （3）、 E ₁₀₀ T の温度で0GPaに達します <サブ>クリティカル =− a / b =3372KこれはTaの融解温度に非常に近い[15]。

a 弾性率と b Taの降伏応力対温度。図5では、弾性率（ E ）の温度依存性のパラメーター化も示しています。 ₁₀₀ ）結果

降伏応力の場合、推奨されるパラメータ化は

ここで Y _ストレス （GPa）と T にあります （ K ）; a =7.99610±0.0415、 b =−0.0039±1.30126×10 ^{− 4} 、および c =7.97381×10 ^{− 7} ±8.32307×10 ^{− 8} 。式から（4）、図5bの実線で示されているように、降伏応力は2次多項式モデルに一致する温度とともに減少する可能性が高いことがわかります。

圧力への依存

「はじめに」のセクションで述べたように、高圧条件下でのTaの熱弾性特性に関して広範な理論的および実験的努力がなされてきました。静的法とは異なり、 E の弾性率の圧力依存性を調べるために、<100>引張荷重による動的法を採用しています。 ₁₀₀ 異なる静水圧下で。

図6は、 E の弾性率の曲線を示しています。 ₁₀₀ 1〜1500Kのさまざまな温度で最大140GPaの圧力に対して。さまざまな温度で動的な方法でさまざまな色のすべての実線が得られますが、正方形のマーカーが付いた点線は、式（1）を使用して静的な方法で得られます。 （1） C の値に基づく ₁₁ および C ₁₂ 0Kで。1Kでの曲線、つまり円のマーカーが付いた赤い実線は、静的な方法で得られた曲線（点線）とほぼ重なっていることは明らかです。これは、現在採用されている動的な方法を示しています。作業は最大140GPaの高圧で適用できます。

さまざまな温度と圧力でのTaの<100>配向の弾性率。これらの曲線の弾性率は E であることを指摘しておく価値があります。 ₁₀₀

図6に示すように、 E の弾性率 ₁₀₀ 600 K以下の温度では、圧力が20から60 GPaに増加するにつれて、下向きの凹面部分が表示されます。 Ruestes etal。 [18]は、MDシミュレーションを使用して、最大60 GPaの圧力の関数として弾性定数を報告しました。結果は、計算された C とよく一致しています。 ₁₁ および C ₁₂ 現在の仕事で。その見返りとして、 E の計算された弾性率 ₁₀₀ また、私たちの結果でも同じ傾向を示しています。しかし、 E の計算された弾性率 ₁₀₀ DFT計算から[33]は圧力の増加とともに徐々に増加し、曲線には下向きの凹面部分は見られません。この不整合が発生する原因は何ですか？一般に、MDシミュレーションで使用されるポテンシャルは、DFT計算と実験結果をフィッティングすることによって得られます。この意味で、DFT計算はMD法よりも精度が高くなります。 Ravelo-EAMポテンシャル[21]は、Ta単結晶のフィッティングDFT状態方程式（EOS）曲線に高圧特性を導入することによって構築されます。フィッティング手順中、EOSのコールドカーブは、格子定数の高次（立方および四次）項を含むように拡張されます。これにより、MDEOSは格子定数の高次項に非常に敏感になります。つまり、この不一致は、Ravelo-EAMポテンシャルが20〜60GPaの圧力下でTaのEOSを正確に記述できなかったことが原因である可能性があります。全体として、弾性率対圧力の曲線は、異なる温度で同様の傾向を示し、弾性率は、圧力が約40GPaを超えると徐々に増加すると結論付けることができます。

結論

この論文では、MDシミュレーションを実行して、<100>引張荷重によるTa単結晶の弾性特性の温度および圧力依存性を調査しました。最初に、0 KでのTaの弾性特性とさまざまな静水圧の観点から、Zhou-EAMとRavelo-EAMを含む2種類のEAMポテンシャルについて比較研究を行いました。結果は、Ravelo-EAMポテンシャルが、さまざまな圧力下でZhou-EAMポテンシャルよりも良好に動作することを示しています。次に、Ravelo-EAMポテンシャルに基づいてTa単結晶の引張挙動のMDシミュレーションを実行します。観察結果は、Taが<100>引張荷重下で破壊前にBCC-FCC相転移を経験し、モデルサイズとひずみ速度がTa単結晶の引張挙動に明らかな影響を与えないことを示しています。さらに、 E の弾性率 ₁₀₀ は、温度が1から1500 Kに上昇すると〜136から〜79 GPaに直線的に減少し、降伏応力は温度が上昇すると〜8から〜4 GPaに減少し、2次多項式に準拠します。最後に、弾性特性の圧力依存性は0〜140 GPaで実行され、観察結果は、全体的な圧力の増加とともに弾性率が増加することを示しています。 MDシミュレーションの結果は、Ravelo-EAMポテンシャルがより高い圧力と E の計算式で適切に動作することも示しています。 ₁₀₀ C を使用する ₁₁ および C ₁₂ 140GPa未満の圧力で。

変更履歴

略語

BCC：

体心立方

DAC：

ダイヤモンドアンビルセル

DFT：

密度汎関数理論

EAM：

埋め込み原子法

EFS：

拡張フィニス-シンクレア

EOS：

状態方程式

FCC：

面心立方

LAMMPS：

大規模な原子/分子の大規模並列シミュレータ

MD：

分子動力学

NPT：

等温定圧集団

NVT：

正準集団

OVITO：

視覚化ツールを開く

PBC：

周期境界条件

Ta：

タンタル