令和３年度高性能汎用計算機高度利用事業「富岳」成果創出加速プログラム「省エネルギー次世代半導体デバイス開発のための量子論マルチシミュレーション」
成果報告書

令和４年５月３０日
名古屋大学・未来材料システム研究所・特任教授
押山

淳

目次

1．補助事業の目的 ・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・2
2．令和 3 年度（報告年度）の実施内容
2-1．当該年度（令和 3 年度）事業実施計画 ・・・・・・・・・・・・・・・・・・・・・・・・・・22-2. 実施内容（成果）(1) 量子論第一原理計算による表面・界面原子機構解明と実験との比較

・・・・・・・・・・・・・3① GaN デバイスに最適なゲート絶縁体の探索

・・・・・・・・・・・・・・・・・・・・・・ 3② GaN/InN 極端ミスマッチ界面での新構造発見

・・・・・・・・・・・・・・・・・・・・・5③ 第 一 原 理 計 算 に よ る G aN 中 の 転 位 と 不 純 物 の 複 合 体 の 構 造 的 同 定・・・・5④ 第一原理計算と機械学習による GaN 成長表面のナノメートルスケール構造の決定・・・・・6⑤ 第一原理計算による SiC 表面ステップ構造の決定と水素吸着エネルギーのエナージェティクス

・・・・・・・・8⑥ 界面欠陥による電子状態高精度計算法開発とその SiC/SiO2 系への応用 ・・・・・・・・・・10⑦ SiC/SiO2 界面における NO アニール後の界面原子構造の解明

・・・・・・・・・・・・・・11⑧ 第一原理分子動力学法による SiC/SiO2 界面の欠陥と NO ガスの反応過程解析 ・・・・・・・13⑨ 3D-NAND フラッシュメモリのメモリ機能発現機構の解明 ・・・・・・・・・・・・・・・・14(2) 統合デバイスシミュレータによる省エネルギーMOSFET の特性解明と実デバイスとの比較 ・・15
(3) エピタキシャル成長統合シミュレーション技術の開発と有用性の検証 ・・・・・・・・・・・・ 18① SiC(0001)成長表面での吸着原子の拡散機構の解明 ・・・・・・・・・・・・・・・・・・・ 18② 窒化物半導体の MOVPE 成長の量子論的シミュレーション ・・・・・・・・・・・・・・・・20③ GaN MOVPE における気相反応のデータ同化 ・・・・・・・・・・・・・・・・・・・・・・22(4) GaN エピタキシャル成長素過程の解明：TOF-MS を用いた MOVPE 気相反応解析 ・・・・・・・・・23(5) プロジェクトの総合的推進 ・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・27(6) 新たな密度汎関数法(DFT)計算手法の開発 ・・・・・・・・・・・・・・・・・・・・・・・・・28① ニューラルネットワークを活用した新たな DFT 計算手法の開発 ・・・・・・・・・・・・・28② 量子コンピューティングでのハイブリッドアルゴリズムの開発 ・・・・・・・・・・・・・30
2-3. 活動（研究会の活動等） ・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・32
2-4. 実施体制 ・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・32
別添１ 学会発表等実績 ・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・34-1-

補助事業の名称「富岳」成果創出加速プログラム
省エネルギー次世代半導体デバイス開発のための量子論マルチシミュレーション

１．補助事業の目的
量子論第一原理物質計算により、アモルファスゲート絶縁膜／半導体に代表されるデバイス界面の特性を解明し、さらに第一原理計算アプリケーションと半導体デバイスシミュレータおよびエピタキシ
ャル成長プロセスシミュレータとを統合した新たなシミュレーション技術を確立する。シミュレーシ
ョン結果と実験結果の検証により、次世代省エネルギー半導体デバイスにおける科学と製造技術の発展に資する。

２．令和３年度（報告年度）の実施内容
２－１．当該年度（令和３年度）の事業実施計画(1) 量子論第一原理計算による表面・界面原子機構解明と実験との比較検討
令和 2 年度においては、GaN パワーMOSFET のゲート絶縁膜材料として有望な、アルミナシリカ・ア
モルファス膜（AlSiO）を MD（Molecular Dynamics）シミュレーションで作成し、本事業参加の名古屋大学、豊田中研実験グループが達成した良好なデバイス性能の物理的起源を明らかにした。また
GaN 中の転位とアクセプター原子結合体の構造と電子特性を明らかにし、実験グループにおける pn接合の漏れ電流の起源を明らかにした。令和 3 年度においては、上記 AlSiO 膜と GaN 界面でのキャ
リヤー・トラップ減少の起源解明、転位と他の不純物原子の複合体の特性解明、を行う。さらに新たな展開として、TEM 実験との共同による GaN/InN ミスマッチ界面構造の決定とその物理的起源の解明、を行う。また SiC エピタキシャル成長の舞台となる、SiC ステップ表面の同定を行う。
(2) 統合デバイスシミュレータによる省エネルギーMOSFET の特性解明と実デバイスとの比較
令和 2 年度は、非平衡グリーン関数法に基づくデバイスシミュレーションにより、本事業参加の東北大学実験グループによる Si ナノワイヤ型 Body-Channel-FET 実デバイスの電流電圧特性を、パラメータなしに再現することに成功した。しかしながら、オン状態での電流値は実験値を超えており、
これは散乱過程の十分な取り込みができていないためだと判明した。令和 3 年度は、キャリヤー・フ
ォノン散乱、キャリヤー・ラフネス散乱の定量的取り込みを行い、より高精度なデバイスシミュレー
ションを行う。
(3) エピタキシャル成長統合シミュレーション技術の開発と有用性の検証
令和 2 年度においては、GaN エピタキシャル成長の素過程を量子論計算で明らかにした。成長表面はGa リッチであり、そこでの Ga-Ga の弱いボンドが成長のホットスポットとなり、窒素が取り込まれていくこと、成長温度では、その表面は 2 次元液体化していること等、驚くべき知見が得られた。令和 3 年度は富岳の本格的利用が見込まれるので、より大規模長時間な MD(Molecular Dynamics)シミュレーションにより、GaN エピタキシャル成長の現実ガス圧、成長温度での機構解明を行う。また GaN
基板成長法として近年注目を集めているアモノサーマル法における素過程解明を開始する。MD シミュレーションにより、最適な成長面方位、鉱化剤選択の処方箋を与えることを目指す。また、SiC のエピタキシャル成長における原子反応過程の探索を行う。

-2-

(4) GaN エピタキシャル成長素過程の解明：TOF-MS を用いた MOVPE 気相反応解析
令和 2 年度においては、精密 TOF（Time-Of-Flight）実験により、GaN 成長炉内でのアンモニア関連の分子状態を明らかにし、その結果、成長炉ガス中ではアンモニアは分解せず、成長表面に到達することを明らかにし、研究項目(3)の成長素過程解明量子論計算につなげた。令和 3 年度においては、
Ga ソース関連分子［TMGa (Tri-Methyl Gallium)およびその分解物］の成長炉内での分子状態解明を行い、ガス中反応の量子論計算と比較検討を行う。
(5) プロジェクトの総合的推進
本事業を効率的かつ効果的に推進するため、プロジェクト推進会議を適宜開催し、参画各機関の連携・調整を実施する。とくに、シミュレーションデータと実験データとの比較、両者のフィードバッ
クによるシミュレーション技術の改善促進、有用な実験データの取得を促進させる。また、プロジェ
クト全体の管理・調整業務を行う。
２－２．実施内容（成果）(1) 量子論第一原理計算による表面・界面原子機構解明と実験との比較検討① GaN デバイスに最適なゲート絶縁体の探索[名大計算グループ + 名大実験グループ + 阪大実験グループ + 富士電機実験グループ]
GaN パワーデバイスのゲート絶縁膜としては、SiO2 と Al2O3 が主に用いられているが、一長一短がある。すなわち、より広いバンドギャップ（9.0 eV）をもつ SiO2 は GaN 界面で電子、正孔ともに十分なバンドオフセットをもつが、誘電率がそれほど大きくない（ ε = 3.8 ）。一方 Al2O3 は大きな誘電率、十分なバンドオフセットが取れない場合（ ε = 8.5 ）を持つが、バンドギャップが狭いため（6.7 eV）がある。名大／豊田中研の実験グループでは、両者の混合酸化物 (Al2O3)1-x(SiO2)x（以下 AlSiO と略記）を用いて GaN(0001)との界面を形成し、良好なデバイス特性を得ている。名大計算グループでは、
AlSiO／GaN 界面の優位性を密度汎関数理論（Density Functional Theory: DFT）に基づく第一原理
計算で調べ、昨年度には界面でのダングリングボンドに起因するキャリヤー・トラップが消失することを発見した。今年度はさらにその消失の成因を明らかにすると同時に、界面特有の遷移相の存在を見出した。
界面構造の決定は、ボルンオッペ
ン ハ イ マ ー (BO) 分 子 動 力 学 法（Molecular Dynamics: MD）による melt-quench スキームにより行われた。界面を図 1(a)に示すスーパ
ーセル構造でモデル化し、このター
ゲットの温度を 6000K まで上昇させ、その温度で 10 ピコ秒(ps)間熱平衡化させた。その後 100K/ps あ
るいは 200K/ps の速度でクエンチすることによってアモルファス
AlSiO を作成した[図 1(b)]。アモル

図１：(Al2O3)1-x(SiO2)x/GaN 界面のスーパーセルモデル。(a) 2 種の結晶から成る初期状態。(b)melt-quench 法によって作成された
アモルファス AlSiO/GaN 界面。

-3-

ファス作成前には界面に無数に存在した Ga ないしはN のダングリングボンドはアモルファス作成の過程で消失し、最終安定構造では、Ga-face の界面では Ga-Oのボンド、N-face の界面では N-Al および N-Si のボンドが形成されている。これはアモルファス作成中に
AlSiO 領域から O、Al、Si 原子が界面に拡散し、ダングリングボンドを終端するためである。実際この界面構造での各原子種の空間分布とその原子種の電荷を図2 に示す。Ga-face 界面に拡散してきた酸素原子はより負に帯電し、また強固な Ga-O ボンドを形成している
こと、N-face 界面に拡散してきた Al、Si は N に電子を供与して強いボンドを形成していることが見て取れる。これは一般的には electron-counting rule という

図 2：AlSiO/GaN 界面での各原子の分布とその原子が有する電荷（電子過多がマイナス
符号）。

化合物系に特有な構造安定化メカニズムである。
得られた安定界面構造での局所状態密度（Local Density of States:
LDOS）を図 3 に示す。界面ダング
リングボンドの消失に対応して、
ギャップ付近には局在した電子状態は存在しない、すなわちキャリヤ
ー・トラップが存在しないことが示されている。
他の外因性（不純物等）
起源のトラップの存在は調べられていないが、
AlSiO の優位性の一端が示された。また図 3 の界面付近の LDOS を詳しく見ると、界面でのバンドオフセットは急峻ではな
く、徐々に変化していくことがわかる。これは Si/SiO2 系でも知られた

図 3：DFT の melt-quench 法で作成されたアモルファス(Al2O3)1x(SiO2)x/GaN 界面構造（上図）と対応する位置での局所状態密度(LDOS)。LDOS の縦軸の原点は GaN の価電子帯上端のエネルギー値であり、白い領域がエネルギーギャップ。

ことであり、半導体／絶縁体界面に特徴的なことである。ナノ構造体の物性解明にとってはひとつの重要なポイントであろう。以上の結果は Journal of Applied Physics 誌への投稿準備中である。
アモルファス SiO2 も重要なゲート酸化膜の候補である。名大計算グループは各大学および富士電機の実験グループとの共同で、アモルファス SiO2/GaN(0001)界面の特性を解明する DFT 計算を開始した。この界面ではしばしばホール・トラップが出現することが知られているが、その起源を解明すること、とくにゲート酸化膜堆積時に不可避的に生成される GaOx 遷移層の効果を明らかにし、富士電機における今後のデバイス開発への指針を与えることが目的である。界面付近の原子構造を melt-quench
法で明らかにしたのが図 4 であり、Ga-O-Si の構造ユニットが出現することがわかった。さらにこの-4-構造では酸素のローンペアに基因するホール・トラップが出現することがわかった。実験結果との突き合わせを行い、SiO2 系のゲート酸化膜形
成の今後の方針を策定した。
② GaN/InN 極端ミスマッチ界面での新構造発見[名古屋大学実験グループ + 名古屋大学計算グループ]
格子不整合のヘテロ構造界面では、生成された膜の厚さが、ある臨界膜厚を超えると界面起因のストレスを緩和するためのミスフィット転位が生じることはよく知られている。しかし格子不整合が極めて大きい場合、界面構造がどうなるかはわかっていない。名古屋大学未来材料・システム

図 4：SiO2GaOx /GaN
(0001)界面での構造とホール・トラップ電子
状態（黄色の雲）。Ga,
Si, O は緑、青、赤球。

研究所五十嵐グループでは、10%の格子不整合が存在する GaN/InN ヘテロ構造で、図 5(a)にあるように、ミスフィット転位の存在しない平坦な界面が出現することを発見した。平面電子顕微鏡写真では、10%の周期の違いによるモアレ縞も観測されている。原子スケールでどのような構造になっているのか、極めて興味深い。図 5(b)は RSDFT コードを用いた DFT 計算による構造探索の結果、得られた安定構造である。平面図ではモアレ状のパターンも得られた。断面透過電子
顕微鏡写真（省略）と計算で得られた断面イメージの比較も良好である。得られた成果は論文準備中である。
図 5：極端格子不整合系 InN/GaN の平坦界面。(a;左) 電顕平面
図。(b;右) DFT による安定構造。平面図（上）と側面図（下）。
赤丸: In、緑丸：Ga、小丸：N。

(b)

(a)
③

第一原理計算による GaN 中の転位と不純物の複合体の構造的同定[名古屋大学計算グループ]

GaN 中には線欠陥（転位）が少なくとも 106cm-2 程度存在し、デバイスの性能に大きな影響を与えている。名大計算グループでは令和 2 年度に、GaN 中のらせん転位の構造を明らかにし、そのらせん転位が Mg アクセプターを引き付け、転位芯の周りの領域で p 型領域が n 型領域に変化してしまうことを見出した。これは名大実験グループにおける pn ダイオードの絶縁破壊現象を見事に説明した。
GaN 中にユビキタスに存在する水素不純物原子は、結晶中では Mg アクセプターと複合体を形成し、Mg アクセプターを不活化することが古くから知られてきた。この不活化には複合体の形成によって生じる電子状態の変化が関わっており、これまでによく研究されてきたが、この水素原子が転位近傍で-5-Mg とともに形成する複合体については明らかになっていない。これを明らかにするため、令和 3 年度
は大規模な DFT 計算を行い、らせん転位、Mg 不純物、H 不純物の複合体構造を明らかにした。
図 6 に計算で用いたモデルを示す。
一般的に、格子欠陥を扱う場合には欠陥で生じる歪みを緩和させるために大きなスーパーセルが必要になる。さらに、らせん転位がある場合の歪みは距離が離れても完全には収束しないので、真空層を置き、擬水素で終端することで結晶中にらせん転位がある様子を表現している。このスーパーセルには
およそ 800 個の原子が含まれている。
このような大きな系の計算には並列計

図 6: スーパーセルの模式図。最も外側の黒枠がセルの境界。緑で示した GaN の周囲に真空層を置き、擬水素で終端している。中心にらせん転位(赤十字)を置いている。このスーパーセルは 6×5×4 個の GaN 結晶の単位胞を含んでいる。

算に特化したプログラムが適している。本研究では密度汎関数理論に基づく
Kohn-Sham 方程式を実空間差分法によって解く RSDFT コードを用いた。さ
らに今回の計算は逆格子空間の原点のみを考慮すればよいので、RSDFT コー
ドの実数バージョンを使うことで計算資源を削減した。図 7 は結晶中の水素とMg 複合体構造および転位近傍での複合体構造である。結晶中では図 7(a)のよ

図 7: 結晶中(a)および転位近傍(b)での Mg および H の安定配置。H:青、Mg:オレンジ、Ga:緑、灰:N を表している。赤十字は転位線の位置を示している。(b)の右図は c 軸方向から見たもので、赤点線は(b)の左図の H 近傍の拡大図に対応した領域を示している。

うに Mg と結合した N に H が隣接し、H-N-Mg の配置をとる。Mg がない場合、つまり単結晶 GaN 中でも同様に H-N-Ga の配置をとる。一方、転位がある場合には転位を挟んで Mg の配置とは反対側になり、H は Mg と離れた配置をとる。さらに H-N-Ga ではなく、N-H-Ga の順に並ぶ。この N-Ga の間に H が配置される現象は、Mg の有無によらず、らせん転位の歪みによって生じていると考えられる。線欠陥と点欠陥（不純物原子）の複合体の構造同定を量子論計算で行った世界初の例である。
④ 第一原理計算と機械学習による GaN 成長表面のナノメートルスケール構造の決定[九州大学計算グループ + 名古屋大学計算グループ]
MOVPE 成長中の GaN 表面は、大気圧下で原料ガスやキャリアガスにさらされている。そのため、真空チャンバー内で観察される表面構造とは異なった、雰囲気分子に被覆された表面が出現し得る。一方で、成長表面は、原料分子さらには不純物分子の吸着・脱離、表面拡散、表面反応といった成長素過程の舞台でもある。例えば、N 原料の NH3 は Ga adatom の Ga-Ga 結合間に取り込まれ、NH として表面拡散する機構が、名大計算グループの研究から明らかにされている。このような素過程レベルから
定量的な結晶成長シミュレータをボトムアップに構築するためには、リアリスティックな表面構造の同定が不可欠となる。従来は(0001)-2×2 という小さい表面周期性の仮定の下で第一原理計算が行われ-6-ており、Ga 原料分圧と H2 分圧に応じて、Ga adatom 表面または H 吸着表面の一方が出現するとされていた。しかし、実際にはこれらの中間組成の成長表面が出現していると考えられる。今回、計算系の表面をナノスケールにまで広げ、成長中の表面構造を量子論第一原理計算で解明した。
ナノスケール表面の構造決定の難しさは、スーパーセルのサイズ増大とともに、候補構造の数（吸着サイトにおける吸着原子の配置パターン数）が爆発的に増加することに起因する。今回取り組ん
だ表面組成においては約 49 万通りの候補があり、網羅的に第一原理計算を実施することは不可能である。近年、機械学習を活用した計算科学が注目を集めている。結晶構造予測の分野で歴史が長いのは進化計算に基づくアプローチであるが、計算コストが高い（1 トライアルの計算時間が長い）大規模系の第一原理計算とは相性が悪い。そこで、より少ないトライアル回数でブラックボックス最適化が可能とされるベイズ最適化と、並列計算性能に

図 8：全 170 トライアルのベイズ探索の履歴。赤点は今回の探索で最も安定な構造、破線は従来のモデルから構成可能な安定な構造のエネルギーを比較対象として示している。

優れた量子論第一原理計算（RSDFT)コードを組み合わせた。ここでは、ベイズ機械学習に基づいて
Ga adatom 表面と H 吸着表面の混合エネルギーが低く予測される構造と、その予測が不確かな構造がバランス良く提案されることにより、効率的な探索が可能になっている。本研究は、大規模表面系
における構造探索の先駆的研究事例となるものである。
図 8 にベイズ探索の履歴を示す。横軸はトライ
アル番号であり、ベイズ最適化により提案された

図 9：GaN MOVPE 中に出現し得る安定な表面構造（Ga adatom 表面と H 吸着表面の混合表面）。

表面構造について実施した第一原理計算の通し番号である。縦軸は Ga adatom 表面と H 吸着表面の混合エネルギーであり、最小化したい量である。回帰モデルは逐次的に更新される。つまり、n 回分のトライアル結果を用いて学習された回帰モデルの予測と不確かさに基づいて、トライアル番号 n+1 の表面構造が提案される。初めはランダムな構造からスタートし学習も進んでいないので、混合エネルギーが 2.5eV を超える不安定な構造が提案されている。25 トライアル程で急速に学習が進み、その後の100 トライアル間は、大きく振動しながらもエネルギー減少の傾向を示している。不安定側に振れている構造は、予測が上手くいかなかった構造、または、不確かさに基づき提案された構造に対応すると考えられる。
図 9 に今回の探索で発見された安定構造を示す（図 8 の赤点に対応する）。この構造は Gaad (2×2)と3Ga-H (2×2)のモチーフ（図中の赤タイルと青タイル）で隙間なく覆われている。すなわち、エレクト
ロン・カウンティング則を局所的に（2×2 表面セルの範囲で）厳密に満たしている。ここで、探索の問題設定として、この制約を課していなかったことに注意が必要である。成長表面がこの特徴を有するこ-7-

とは必ずしも自明ではなく、新たな知見と言える。また、青タイル内の H 吸着位置はそれぞれ異なっている。これはタイル間の相互作用を緩和するように配置された結果だと考えられる。図 9 の他に見つかった安定構造も同様の特徴を有していた。今回の探索で計算された構造は全候補構造のごく一部であるので、同様に安定な構造は多数あると考えられるが、それらが有するであろう上述の特徴を見出すことができた。
⑤ 第一原理計算による SiC 表面ステップ構造の決定と水素吸着エネルギーのエナージェティクス
[名古屋大学計算グループ]
SiC のエピタキシャル成長は、SiH4、C3H8 をソースガスとして用いる CVD(Chemical Vapor Deposition)法で行われる。さらに膜質向上と成長速度増進のために HClを加えた HCVD (Halide CVD)法もしばしば使われる。
昨年度の報告書で説明したように、成長炉内気相反応に対する我々の第一原理計算によれば、気相中には SiCl2,
C2H2 および H2 がマジョリティとして存在する。それらが成長表面に到達し、原子反応を経て表面ステップ端から、エピ膜に取り込まれるのがステップフローの成長である。その機構解明のためには表面構造、とくにステップ
端の構造、そこでの水素被覆状態を知る必要がある。今年度、そのステップ構造の決定とエナージェティクスの解明を行った。
図 10 はステップ端構造の原子配置の模式図である。
(0001)面の成長では、<11-20>方向に微傾斜した基板が使

図 10: SiC(0001)面のステップ構造の模式図。中央を上方テラス、周囲を下方テラスとし、5 種類のステップ端、SC, Si2, C2,
Si3,C1 が示されている。Si が青、C が茶、
ステップ端 Si が緑、ステップ端 C が黄で示してある。

われる。その場合には Si と C がともにステップ端に現れる SC ステップ端が出現する。またステップ
がキンク状になる場合には、局所的に<1-100>方向のステップとなるので、Si2、Si3、C1、C2 という4 種類のステップ端構造が可能である（図 10: ステップ端の原子種と構造緩和前のダングリングボン
ドの数で、4 種類をラベルしている）。そのステップ端に水素原子は何個どのように吸着するかの可能性は数十種類あるが、今回その全ての可能性を DFT 計算で網羅的に探索し、その結果 30 種類の安定
ないし準安定吸着構造があることを見出した。そのすべての（準）安定構造での吸着エネルギーが表 1に示してある。C と Si ステップ端を比べると C 端への水素吸着エネルギーの方が大きい。これは SiH ボンドより C-H ボンドの方が強いためである。5 種類の異なるステップ端それぞれの場合に、最も吸着原子 1 個当たりのエネルギーの得が大きい構造は、SC[H2B]、C1[H2]、C2[H2]、Si2[H2A]、
Si3[H2A]であり、その吸着エネルギーはそれぞれ、1.40 eV、2.33 eV、1.46 eV、1.51 eV、0.66 eV となっている。Si3 構造での水素吸着エネルギーが顕著に小さいことが見て取れる。安定吸着構造の例として図 11 に SC[H2B]、Si2[H2A]、Si3[H2A]構造を示す。一般的に上方テラスの Si エッジ原子に水素が吸着した構造が安定であることが見いだされた。
水素吸着エネルギーが求まったので、気相中の水素分子と成長表面が局所熱平衡にあるとすると、各-8-

ステップ端での水素の被覆率を計算することができる。結果が図 12 である。ステップ端のみならず表面テラスでの被覆率も求めた。テラス上での被覆率に比べてステップ端での水素被覆率は低く、ステッ
プ端が成長反応のスポットになることが示された。ただし、Ｃ1 ステップの被覆率は極めて高く、このステップでの反応が律速となることが予想される。以上のステップ構造の同定により、SiC-CVD 成長
反応機構解明の舞台が整った。この成果は、Physical Review B 誌に投稿中である。
表 1：SC、Si2、Si3、C1、C2 各ステップ端に SC、Si2、Si3、C1、C2 各ステップ端に水
0素原子が吸着した場合の 30 種の吸着構造とその吸着エネルギー ∆Ead 。吸着した水素原子
数 nH
0
ad

∆E

で吸着構造がラベルされている。水素原子 1 個当たりの吸着エネルギーはnH 。

Step edge

SC

Si2

Si3

C1
C2

H1-A
H1-B
H1-C
H2-A
H2-B
H2-C
H3-A
H3-B
H4
H1-A
H1-B
H2-A
H2-B
H2-C
H2-D
H3-A
H3-B
H4
H1-A
H1-B
H2-A
H2-B
H2-C
H3-A
H3-B
H4
H1
H2
H1
H2

Number of H adsorbed

𝒏𝒏𝐇𝐇
1
1
1
2
2
2
3
3
4
1
1
2
2
2
2
3
3
4
1
1
2
2
2
3
3
4
1
2
1
2

Adsorption energy [eV]𝚫𝚫𝑬𝑬𝟎𝟎𝐚𝐚𝐚𝐚

0.52
−0.09
1.20
1.23
2.80
1.87
1.76
0.89
3.28
1.40
0.29
3.01
1.90
1.32
0.82
1.88
0.99
1.37
0.40
0.05
1.32
1.01
0.02
0.86
0.14
0.27
2.30
4.65
0.51
2.91

-9-

𝚫𝚫𝑬𝑬𝟎𝟎𝐚𝐚𝐚𝐚 /𝒏𝒏𝐇𝐇
0.52
−0.09
1.20
0.62
1.40
0.94
0.59
0.30
0.82
1.40
0.29
1.51
0.95
0.66
0.41
0.63
0.33
0.34
0.40
0.05
0.66
0.51
0.01
0.29
0.05
0.07
2.30
2.33
0.51
1.46

図 11：SC、Si2、Si3 ステップ端での最も安定な水素吸着構造[側面図（左）と平面図（右）]。
元素を示すカラーコードは図 10 と同じ。

⑥ 界面欠陥による電子状態高精度計算法開発とそのSiC/SiO2 系への応用[東京工業大学計算グループ +Quemix 計算グループ]
半導体デバイス開発には、デバイス中の点欠陥特性を
荷電状態まで含めて正しく理解することが必須である。
本年度は、デバイス中の高精度な点欠陥計算手法の開発に着手した。高精度な第一原理点欠陥計算を実行するには、計算システムサイズに対する有限サイズ補正を適切に行う必要がある。これまでに、半導体結晶中における点欠陥の高精度な電子状態計算手法は既に確立したもの

図 12：各安定ステップ端構造での水素被覆率の温度依存性。摂氏 1625 度は典型的なエピ成長温度。

が存在する。一方で、本研究課題が着目する表面・界面における点欠陥の有限サイズ補正法は今だに確立したものが存在しない。本年度は、表面界面における点欠陥の有限サイズ補正法の開発に取り組んだ。第一原理計算では、周期的境界条件を用いた計算が行われ、帯電した点欠陥が周期的に整列することになる。本来求めたい、孤立極限の点欠陥の形成エネルギーを求めるためには、それら整列した点欠陥同士の相互作用を分離し、消去・補正を行わなくてはならない。今回、周期的に並んだ表面
点欠陥のスラブモデルに対して、ポアソン方程式を解き、点欠陥同士の静電相互作用力を見積もることから、点欠陥間の相互作用補正を実施した。開発した表面点欠陥補正法を、これまで理論・実験的に同定されている SiC/SiO2 界面の PbC センター（炭素 adatom のダングリングボンド）に対して適用を行った。その結果、(0/-)（荷電状態が 0（中性）と-1 価とで、エネルギーの安定性が切り替わるエネ
ルギー準位）が価電子帯から 2.5 eV（速報値）高いエネルギー位値となった。この計算結果を元に、実験グループとの比較検討を行なっている。この成果により、表面・界面における点欠陥の高精度電子状態計算手法が確立されるものと考える。
また、昨年度・本年度に亘って行なってきた研究課題の 1 つとして、SiC/SiO2 界面に形成された単一光子光源（点欠陥の一種）に対する電子状態解析を行なった。SiC/SiO2 界面に形成された点欠陥の中には単一光子光源として振る舞うことが報告されているものがある。また、室温において ODMR
（光検出磁気共鳴）測定が可能であり、またその高い輝度から量子センサとして大きな注目を集めている。しかし、SiC/SiO2 界面に形成された単一光子光源（以下、表面 SPS と呼ぶ）は、表面の位置によってその発光波長にばらつきがあることが報告されており、そのばらつきの発生メカニズムは未解明であった。我々は、第一原理計算を用いて、その発光波長のばらつきメカニズムを解明した。計算- 10 -

としては、第一原理分子動力学計算（MD）により表面 SPS サンプルを 5 つ用意し、その発光波長ばらつきを観察した。理論計算の結果、ばらつきとしては 0.3 eV 程度のばらつきを見出すことができ、実験においても、同程度の発光波長ばらつきが報告され、定性的にも実験と矛盾しない計算結果を得ることができた。また、ばらつきの微視的なメカニズムが局所的な歪みであることを見出した。得られた結果を実験へとフィードバックし、実際に以下のような実験を行い、計算結果の検証を行った。
実験で得られた表面 SPS に対して、アニールを行い再度同じ表面 SPS を観測することにより、発光波長を再度観察した。アニールによる局所歪みの変化が期待され、それに応じて発光波長のズレを観察されるものと期待した。実際に、実験の結果アニール後に発光波長にズレが生じていることが確認され、理論計算との整合性が確認された。得られた結果は論文として出版を行なった。[Applied
Physics Letters 118, 204005 (2021)]⑦ SiC/SiO2 界面における NO アニール後の界面原子構造の解明[神戸大学計算グループ ＋

北海道大学計算グループ]

SiC は次世代パワーデバイス用チャネル材料として期待されているものの、SiC-MOSFET の移動度は、SiC
バルクの移動度に比べて極めて低く、材料の性能を活かしきれていない。酸化プロセス後の NO ガスアニー
ルにより界面欠陥密度を低減することで、移動度が向上するという報告がある。しかし、アニールガス原子の移動度向上に対する役割やアニール後の界面原子構造は明らかになっていない。2020 年度は、計算量縮減のため酸化膜を含まないバルクモデルで NO ガスアニール後の窒化物の原子構造・電子状態を調べ、a 面
が NO ガスアニールに最も適した MOS 界面であるという知見を得た。2021 年度は、SiC/SiO2 の酸化膜の影
響を考慮するため、バルクモデルよりも原子数が多い界面モデルを用いて窒化物の原子構造・電子状態を解析した。a 面、m 面、Si 面、C 面の 4 種類の SiC/SiO2 界面に対し、NO ガスアニール前後の原子構造と全エネルギーを比較した。SiC の C 原子が SiO2 の O 原子と CO 結合する構造も存在しうるため、a 面、m 面、C面界面では、CO 結合を有する界面と有さない界面も比較した。図 13(a)に a 面界面の計算に用いた CO 結合を有さないアニール前の界面原子構造を、図 13(b)にアニールにより SiO2 直下の SiC 層(top)に窒化膜が形成された界面原子構造を示す。アニール温度は 1000K、NO ガスの分圧は 1atm、反応複製生物である N2ガス分子と CO ガス分子は、それぞれ 0.25×10−1～0.25×10−5atm、10−1～10−5atm の範囲で化学ポテンシャ
ルを変化させ、窒化膜の形成
エネルギーを調べた。
図 14 に窒化物の形成エネ
ルギーを示す。a 面を用いた
SiC/SiO2 界面に形成される窒化膜は、m 面、Si 面、C 面の界面に形成される窒化膜よりも形成エネルギーが低い。このように、界面モデルを用いた計算でも、バルクモデルを用いた計算結果と同様に、a 面

図 13: 窒化膜形成計算モデル。(a)アニール前。(b)アニール後に第 1 層に窒化膜形成。(c)アニール後に第 2 層に窒化膜形成。青球は Si 原子、茶球は C 原子、赤球は O 原子、緑球は N 原子、灰球は H 原子である。図は a面界面の例である。

- 11 -

図 14: CO 分子と N2 分子の分圧に対する窒化膜形成エネルギーの変化。(a)界面に CO 結合を有さない場合。(b) 界面に CO 結合を有する場合。形成エネルギーの低い方が窒化膜を形成しやすい。

を用いた MOS 界面が NO ガスアニールに適しているという結論を得た。次に、SiC/SiO2 界面に原子レベル
で平坦な窒化膜層が形成されることを確認するため、窒化膜が SiC/SiO2 界面の SiO2 直下の SiC 第 1 層(top)に形成される場合と、図 13(c)のように界面遷移層を挟んだ SiC 第 2 層(2nd)に形成される場合の形成エネル
ギーを比較した。その結果、図 14 に示すように、CO 結合を有さない界面では top 層の形成エネルギーの方が低く、窒化膜は SiC/SiO2 界面に局在することが分かった。これは、NO ガスアニールにより導入される N 原子が界面に局在するという実験結果と符合する。一方、CO 結合を有する界面では、遷移層を介して第 2 層に窒化膜が形成されやすい。第 1 層に窒化膜を形成する構造では、N と O の結合が界面に作られるが、Si、C、N、O の中で電気陰性度の高い N と O の結合は電気的に反発しあい不安定となるため、形成エネルギーが高くなる。しかしながら、遷移層を持つ NO アニール後界面の報告例は、我々の知る限りない。そこで、遷移層を介して第 2 層に窒化膜が形成された構造と図 13(c)のように CO 結合を有さない界面で遷移層を挟まない構造の形成エネルギーを比較したとこ
ろ、遷移層は容易に酸化され、遷移層を挟まない構造に変化することが分かった。したがって、CO 結合を有する SiC/SiO2 界面に形成される遷移層は、アニール時に酸化され消滅するものと予想される。
最後に、図 15 に界面及び基板深部の局所状態密度を示す。SiC/SiO2 界面における CO
結合は、SiC の禁制帯中に欠陥準位を作る。
界面における CO 結合が移動度に及ぼす影響は未解明であるが、NO ガスアニールは、このような禁制帯中に準位を作る構造を滅する。以上のように、NO ガスアニールにより、キャリヤー・トラップ密度低減を期待できる原子レベルで平坦な窒化膜が SiC/SiO2 界面に形成されることが明らかになった。

図 15: 界 面 及 び 基 板 深 部 の 局 所 状 態 密 度 。 太 線 は
SiC/SiO2 界面もしくは SiC/窒化膜界面の SiC 側。細線は、基板深部。CO 結合を有する a 面界面では、SiC の禁制帯
中に準位を作る。

- 12 -

⑧ 第一原理分子動力学法による SiC/SiO2 界面の欠陥と NO ガスの反応過程解析（富士電機計算グループ）
パワー半導体と絶縁膜の界面および絶縁膜中に含まれる欠陥の抑制・制御は、キャリア移動度や信頼性の向上に重要である。しかし絶縁膜まで含む界面の計算モデルは多くの原子数を必要とし、一般的な第一原理計算では計算負荷が大きいため、プロセスガスとの反応の計算科学的な分析が困難であった。今年度は SiC/SiO2 界面において NO ガスによる熱処理による欠陥抑制メカニズムを大規模第一原理分子動力学法により検討した。
東京大学で開発された実空間密度汎関数法により、SiC/SiO2 界面に NO ガスを挿入した系において、第一原理分子動力学計算を行い、界面に存在する欠陥と NO ガスの反応を模擬し、欠陥の構造変化を解析した。SiC は界面の終端構造がよく知られた C 面（000-1）と Si 面（0001）の 2 つを用いた。初めに結晶 SiC と結晶 SiO2 の界面を含む真空スラブモデルを用意して、界面から離れた部分を固定しながらメルトクエンチ法により界面部分を緩和した。これにより界面で SiC 側表面の未結合手
や C=C 結合など、SiC/SiO2 界面において重要と考えられている欠陥を含む界面構造が得られた［図16(a)］。このような構造について、NO ガスを欠陥の近くへ配置して 1600 K 1ps の第一原理分子動力学計算を繰り返した。
C 面界面においては、NO ガスが SiC 側表面に存在した C を SiO2 側へ押し上げて N が C の下へ潜。こり込み、N が SiC 側の Si と 3 つ結合をもつ Si3N 構造がすぐさま（2ps）形成された［図 16(b)］のような N による C 面表面の終端は、SiC 表面の未結合手を低減し、移動度を向上させる一因として考えられており、実験的にも観測されている。これらのことから本計算が SiC/SiO2 界面の欠陥と NO
ガスの反応を定性的に再現することが確かめられた。
C 面において Si3N 終端構造は形成されたが、界面に存在する C=C 結合は NO ガスを数個追加しても分解されず、NO は N と O に解離し、N が C と Si の結合の間に入るような欠陥が生成した［図16(b)の O,C,N からなる欠陥］。この結果から、C 面において NO は C=C の結合を直接的には壊さないことが示唆された。
Si 面界面おいても、界面で Si の未結合手を低減させるような Si3N 構造が形成されたが、周囲に Cが存在しない場所に限られた［図 17(a)］。すなわち、C 面で見られたような N が C を押しのけて Si3Nを作るような挙動は見られず、N は Si に加え C とも結合をもち続け、Si3N 結合が形成されにくい傾向にあることが示唆された。これは、結合相手の Si の位置が C 面では全て SiC 側にあり、N の収まるべき位置が決まっている一方、Si 面では結合相手の 2 つの Si が SiO2 側に存在し、構造の自由度が大きいため、C との結合を切りにくいことが影響していると推察される。これらのことは、C 面のほうが Si 面よりも N の取り込み量が多いという実験結果とも対応する結果であり、面方位依存性を定性的に再現する計算結果であることが確かめられた。
最後に、Si 面において SiO2 側に存在した C=C 欠陥の付近に NO ガスを配置し、その生成物を調べた［図 17(b)］。NO ガスを 1ps 毎に 1 個ずつ追加していくと、C=C 欠陥は C=C 結合を維持したまま、
C-Si の結合間に N が入り込み、Si-N-C-C-N-Si のような結合が形成された。さらに NO を追加すると、C=C の結合間にも N が入り込み C=C 結合は解消された［図 17(b)］。この間、CO や CO2 として脱離する C は見られなかった。これは、C 面でも見られたように NO ガスが C=C 結合を直接的には- 13 -

壊さないことを示唆する結果であり、NO ガスによる C=C 欠陥の分解は Si-N-C の結合形成から始まる多段階の反応である可能性があることが分かった。

図 16：(a) SiC(C 面)/SiO2 真空スラブモデル。
(b)NO ガス導入と 1600K1ps の動力学計算を 2回実施後に形成された Si3N。

図 17：Si 面のスラブモデルにおいて(a)Si3N が形成される場合と(b)C=C 欠陥と NO を反応させた場合の構造。黄、灰、赤、青が Si, C, O, N
原子を表している。

⑨ 3D-NAND フラッシュメモリのメモリ機能発現機構の解明（名古屋大学計算グループ ＋ 東京工業大学計算グループ）
3D-NAND フラッシュメモリはデータ爆発に対応する
不揮発性メモリの大容量化にとって非常に重要であり、省エネルギー半導体デバイスで第一義的な位置を占め、その市場規模は１０兆円に及ぶ。このメモリは電子をSiN 膜中の欠陥に注入することが書込み動作であり、電子を SiN 膜から引き抜くことが消去操作に対応する。しかし、メモリ機能を担っている SiN 膜中の欠陥の原子レベルの構造については明確な理解がなされていない。実験的には Si リッチの SiN 膜を形成することで良好なメ
モリ動作をすることがわかっているので、SiN 膜中の N

H

1.67eVH

図 18: 中性状態の N 空孔と H が分離しているときと複合体を形成したときの原子構造とエネルギー利得の模式図。

空孔関連の欠陥がメモリ機能を担っているという予想はたてられている。2021 年度は SiN 中の N 空孔と N 空孔-H 複合体についてその構造安定性と電子構造を研究した。その結果、図 18 に示すように中性状態では N 空孔と H が分離しているよりも N 空
孔-H 複合体を形成した方が 1.67eV エネルギー的に得をすることを明らかにした。さらに、SiN 中に存在する H 原子は系の荷電状態に依らず N 空孔と複合体を形成した方がエネルギー的に安定となることも確認した。このように、3D-NAND フラッシュメモリの機能を担っているのが N 空孔-H 複合体であることを世界に先駆けて明らかにした。次年度はこの知見をもとに、アモルファスの SiN 膜に- 14 -

ついても研究を行う。
(2) 統合デバイスシミュレータによる省エネルギーMOSFET の特性解明と実デバイスとの比較[大阪大学計算グループ]
非平衡グリーン関数法（NEGF 法）に基づく量子論デバイスシミュレータを用いて、省エネルギー型ナノワイヤトランジスタやナノシートトランジスタなどのナノデバイスの開発指針を早期に得るためのデバイスシミュレーション環境の構築を進めている。本年度は、変形ポテンシャル相互作用ハミル
トニアンを用いたフォノン散乱過程、および、面内分極電荷の効果を量子輸送シミュレータに導入した。
東北大学の遠藤教授らが開発した直径 60 nm、チャネル長 100 nm の縦型ナノワイヤトランジスタ[1]をシミュレーション可能な大規模量子輸送デバイスシミュレータの開発を行っている。有効質量近似を用いて電子状態を記述し、輸送部分に効果的な並列計算手法を用いた大規模 NEGF デバイスシミ
ュレータの開発を行っている。はじめに、フォノン散乱などの散乱過程を無視した弾道輸送条件のもとでシミュレーションを行い、実デバイス特性を概ね良好に再現することを確認した。しかし、オン電流を大きく見積もる問題があった。そこで、この不一致を解消することを目指し、昨年度は、単一の光学
フォノンモードを仮定したフォノン散乱モデルを導入した。この簡単化したフォノン散乱モデルには、
フォノン散乱強度を表すフィッティングパラメータ𝐵𝐵𝜈𝜈 が必要であり、その値を、低電界電子移動度がバルクシリコンの実測値と一致するという条件から決定した。その結果、フォノン散乱を導入すると、オン電流が減少し、その減少割合は、低いソース・ドレイン間電圧の方が大きいことなどがわかった。
しかし、シミュレーション結果は、依然として、測定結果を大きく見積もるという問題が残った。そこで、本年度は、フィッティングパラメータを必要としない、より詳細化したフォノン散乱モデルを導入した。さらに、面内分極電荷効果の導入を行った。
フォノン散乱モデルとして、音響モード変形ポテンシャル散乱（弾性散乱）と、バレー間フォノン散乱（非弾性散乱）とを導入した。音響モード変形ポテンシャル散乱の相互作用ハミルトニアンは、つぎの形を仮定した。
𝐻𝐻ac = 𝐷𝐷ac div𝜑𝜑; 𝜑𝜑 = �𝒒𝒒𝒒𝒒
�𝑏𝑏𝒒𝒒 𝑒𝑒 𝑖𝑖𝒒𝒒⋅ 𝒓𝒓 + 𝑏𝑏𝒒𝒒+ e−𝑖𝑖𝒒𝒒⋅ 𝒓𝒓 ��2𝜌𝜌𝜌𝜌𝜌𝜌 𝑞𝑞
1

ここで、ρ は密度、𝑢𝑢 は音速である。変形ポテンシャル定数 𝐷𝐷ac の値は、MOS 構造 [2] およびナノワイヤ構造 [3] の解析で用いられた文献値(14.6 eV)を用いた。一方、バレー間フォノン散乱の相互作
用ハミルトニアンは、つぎの形を仮定した。
𝐻𝐻𝑗𝑗 = 𝐷𝐷𝑗𝑗 Φ𝑗𝑗 ; Φ𝑗𝑗 = �𝒒𝒒1

�2ρω𝑗𝑗+ −𝑖𝑖𝒒𝒒⋅𝒓𝒓�𝑐𝑐𝒒𝒒𝑗𝑗 𝑒𝑒 𝑖𝑖𝒒𝒒⋅𝒓𝒓 + 𝑐𝑐𝒒𝒒𝑗𝑗𝑒𝑒�

3 種類の g タイプ散乱と 3 種類の f タイプ散乱とを考慮し、それぞれの変形ポテンシャル定数 𝐷𝐷𝑗𝑗 およびフォノンエネルギー ℏ𝜔𝜔𝑗𝑗 は、文献 [4] の値を用いた。

図 19 にシリコンナノワイヤにおける移動度の電子密度依存性を示す。昨年度のモデルでは、相互作用強度をバルクシリコンの移動度(1,350 cm2/Vs)を再現するよう調整しており、低電子密度領域で近い値を示すことがわかる。一方、本年度のモデルは、低電子密度領域における電子移動度が 800 cm2/Vs- 15 -

程度に低下することがわかった。これは、本研究で用いた音響モード変形ポテ
ンシャル定数の値が、MOS 構造・ナノワイヤ構造の解析に用いられている値であるためである。
図 20 の赤丸印に，𝑉𝑉SD = 0.05 V における伝達特性を示す。比較のため、測定値[1]を青丸印で、弾道輸送条件における計算結果を黒実線でプロットした。また、𝑉𝑉G = 0.4 V における電流スペクトルの計算結果を図 21 に示す。本年度のモ
デルは、図 19 に示した通り、バルクシ

リコンの電子移動度より低い電子移動度を示すモデルであるため、昨年度のモ

図 19：シリコンナノワイヤにおける電子移動度の電子密度依存性．赤丸印が変形ポテンシャル相互作用を用いた結果、青丸印は，単一のフォノンモードを仮定した結果。

デルと比べて、オン電流が、より測定値に近づくことがわかった。
図 20 に示した通り、測定結果と比較して、フォノン散乱を導入しても、シミュレーション結果の方が測定結果より若干大きいことがわかる。直径 1 nm、チャネル長 10 nm 程度の非常に細い短チャネ
ルナノワイヤの場合、非平衡分極電荷を考慮すると、定数の誘電率を仮定した場合に比べ、ドレイン電流が小さくなることを確認している。そこで、本年度、非平衡分極電荷の影響を取り入れるためのプログラムを開発した。

図 20：直径 60 nm、チャネル長 100 nmのシリコンナノワイヤトランジスタの伝達特性。赤丸印がフォノン散乱を考慮した結果、黒実線は弾道輸送条件の場合。測定結果を青丸印に示す。

図 21：𝑉𝑉G = 0.4 V における電流スペクト
ルの計算結果。赤実線は、フォノン散乱を考慮した結果、黒破線は、弾道輸送条
件の場合。

非常に細い短チャネルナノワイヤの場合、ソース・ドレイン輸送方向のポテンシャル変化が、輸送方向に垂直な面内方向のポテンシャル変化に比べて大きいため、輸送方向の非平衡分極電荷分布をシミ- 16 -

ュレータに取り入れることが重要であった。一方、遠藤教授らが開発したナノワイヤトランジスタの場合、直径が太い
ため、面内方向の非平衡分極電荷分布を取り入れることが重要である。
はじめに、MOS 構造において非平衡
分極電荷分布を計算するプログラムを開発した。Si/SiO2 界面に垂直な 𝑧𝑧 方向の 1 次元電子状態を自己無撞着的に計算するプログラムを作成した。界面から20 nm (214 原子)の範囲における分極電
荷分布 ρv (z) を、RSDFT を用いてつぎの式により計算した。

ρ𝑣𝑣 (𝑧𝑧) = 2 � δ�Ψ𝑛𝑛𝑛𝑛 �

図 22：MOS 構造における分極電荷分布。灰色実線が RSDFT
による結果であり、それを平滑化したものを赤実線に示す。
バルクシリコンの誘電率を仮定した平均場近似（MFA）の結果を青点線でプロットした。

v.band

ここで、𝛿𝛿|Ψ𝑛𝑛𝑛𝑛 | は、閉じ込めポテンシャルによる波動関数変化を表す。一方、反転層電子密度は、有効質量近似を用いて計算した。計算した分極電荷分布を図 22 に示す。灰色の実線が RSDFT による結果である。比較のため、RSDFT の結果を平滑化したものを赤実線で、バルクシリコンの誘電率を仮定
した平均場近似（MFA）の結果を青点線で示す。界面 (𝑧𝑧 = 0) 付近の極めて薄い領域を除いて、RSDFTの結果が、バルクシリコンの誘電率の結果を再現することを確認した。

(a)

(b)

図 23：直径 30 nm、チャネル長 60 nm のシリコンナノワイヤトランジスタのオン状態におけるポテンシャル分布。(a) 面内動径方向分布、(b) 輸送方向分布。赤実線は、RSDFT を用いて非平衡分極電荷分布を考慮した場合、黒破線は、バルクシリコンの誘電率を用いた平均場近似（MFA）の結果。

つぎに、分極電荷分布の効果を取り入れたデバイスシミュレータのプロトタイプを開発した。図 23に、ターゲットデバイスよりサイズが小さい、直径 30 nm、チャネル長 60 nm のナノワイヤトランジ
スタのオン状態におけるポテンシャル分布の計算結果を示す。分極電荷を考慮した場合、Si/SiO2 界面
付近の面内方向ポテンシャル分布にわずかな違いが生じ（図 23(a)）、チャネル中央付近の輸送方向ポ- 17 -

テンシャルが高くなることがわかる（図 23(b)）
。これらの効果により、オン電流が低下する。図 24 に直径 30 nm、チャネル長 60 nm のナノワイヤトランジ
スタの伝達特性のシミュレーション結果を示す。今後は、直径が 60nm、チャネル長が 100 nm のナノワイヤトランジスタにおけるシミュレーションを行えるようにプログラムの改良を行う。
[1] T. Imamoto et al., Jpn. J. Appl. Phys., 54,
04DC11 (2015).
[2] D. Esseni et at., IEEE Trans. Electron Devices,50, 2445 (2003).
[3] S. Jin et al., J. Appl. Phys., 99, 123719 (2006).

図 24：直径 30 nm、チャネル長 60 nm のシリコンナノワイヤトランジスタの伝達特性。黒実線および青丸印は、平均場近似の結果、赤丸印
は、RSDFT を用いた場合。実線が弾道輸送条件、丸印がフォノン散乱を考慮した結果。

[4] C. Jacoboni and L. Reggiani, Rev. Mod. Phys.,55, 65 (1983).
(3) エピタキシャル成長統合シミュレーション技術の開発と有用性の検証① SiC(0001)成長表面での吸着原子の拡散機構の解明 [名古屋大学計算グループ]
エピタキシャル成長、とくにステップフロー・モードでの成長の素過程は、i) 気相から成長平面への反応分子の飛来、ii) 反応分子の表面テラス上での分解反応、iii) 分解生成物のテラス上およびステップ端近傍での拡散、iv) ステップ端からの結晶薄膜への構成元素の取り込み、の４つの過程に分類できる。今年度、SiC の CVD 成長における iii)の拡散を DFT 計算で調べた。実際の CVD 成長における表面ステップ構造は、上述の研究項目(1)⑤で検討中であり、また ii)の分解反応は令和 4 年度の課題であるが、ここでは、すでに我々の DFT 計算で明らかになった水素原子に被覆されていない SiC(0001)
表面上での、究極の分解反応生成物である Si、C、H アド原子（adatom）の拡散を調べた。

図 25：SiC(0001)表面テラス上での Si、C、H アド原子の表面位置に対する全エネルギーの等高線図。エネルギーの単位は eV、エネルギー原点はそれぞれのアド原子の最安定位置での吸着エネル
ギー。安定吸着位置 T4、H3 が示されている。

計算は表面をスラブ模型でシミュレートし、RSDFT コードを用いて行った。拡散経路の決定と対応する活性化エネルギーの導出は、各アド原子の表面(x,y)座標を固定し、他の自由度を構造緩和させ、得- 18 -

られた全エネルギーの等高線図（landscape）により行われた（図 25）。表面テラスの原子構造は 3x4のスーパーセル模型での安定構造計算により、2x1 の表面原子 buckled 構造が得られたので、それを構造ユニットとして用いた。ステップ構造は以前の我々の DFT 計算（Seino, Oshiyama: APEX, 2020;
PRB 2020）の結果を用いた。
図 25 は、各アド原子に対する表面テラス上での全エネルギー等高線図である。安定な吸着位置とし、H3 サイト（表面原子六角形の真中）があることがわかる。
ては、T4 サイト（表面下 C 原子の真上）Si アド原子の最安定位置は T4(2)サイトであり（カッコ内の数字は飛び出た表面 Si 原子とのボンドの数を示している）、得られた拡散経路は、[1-21]（六方晶では[1-100]）方向には・・T4(2)・H3(2)・T4(1)・H3(1)・T4(2)・・、[10-1]（[11-20]）方向には・・T4(2)・H3(1)・T4(2)・H3(1)・T4(2) ・・である。対応する拡散の活性化エネルギーは 0.65 eV ～ 0.68 eV であり、Si アド原子のテラス上での拡散は等方的であることがわかった。C アド原子の吸着位置、拡散経路は Si アド原子と同じであることが判明した。拡散が等方的であることも同様である。しかしそのエネルギーは大いに異なり、拡散の活性化エネ
ルギー1.45 eV と求められた。H アド原子の吸着位置は T1(1)と T1(2)であり、拡散活性化エネルギーは
1.15 eV～1.35 eV と求められた。ここで得られたテラス上での拡散特性は、キネティックスに対する
今後のシミュレーションの土台となるであろう。
ステップ構造については、以前の我々の計算で詳細な構造が明らかになっている。またその生成エネ
ルギーも求められており、Si2 ステップ（図 10 および図 26(a)参照）が最も出現しやすい構造であることがわかっている。そこでこの Si2 ステップ近傍での全エネルギー等高線図を求め、拡散経路と活性化エネルギーを調べた。等高線図の一例として、Si アド原子に対するものを図 26(b)に示す。テラス上
での等高線図（図 25）とはかなり趣を異にしている。

図 26：(a)SiC(0001)表面 Si2 ステップ構造の側面図。青、緑は Si 原子、茶、黄は C 原子を示す。(b) Si2 ステップ近傍（左が上方テラス、右が下方テラス）での Si アド原子の面内位置に対する全エネルギーの等高線図。

図 27 は決定された拡散経路とその経路に沿った全エネルギーの変化である。ステップ端近傍上方テラスでの安定位置としては H3u、T4u、T1u などがあり、下方テラスでは H3d が安定である。これらを出
発点あるいは到達点とする拡散経路が図 27 に示されている。Si アド原子については、ステップ近傍ではテラス上に比べ、顕著な活性化エネルギーの増加が見られる。これは従来の現象論において、EhrlichSchwoebel (ES)効果として導入されてきたものであり、その量子論計算による裏付けとなっている。
しかし、C アド原子にたいしてはそのような ES 効果はほとんど見られないことが初めてわかった。逆にわずかではあるが、活性化エネルギーの減少も見られ、逆 ES 効果とも言える。いずれにしても、今回の計算で判明した拡散経路活性化エネルギーの定量的決定は、成長素過程の解明のプロトタイプと- 19 -

なるものであろう。以上の結果は Applied Surface Science 誌に発表された。

図 27：SiC(0001)表面 Si2 ステップ近傍での各アド原子の拡散経路（下図）とその経路に沿った
全エネルギーのプロファイル（上図）。

② 窒化物半導体の MOVPE 成長の量子論的シミュレーション[名古屋大学計算グループ + 九州大学計算グループ]
2020 年度は GaN の MOVPE 成長の気相反応の原子レベルの過程を明らかにしたが、2021 年度は、GaN と合わせて、次世代通信技術であるポスト 5G-6G 規格において重要な材料となる窒化物半導体
AlN の MOVPE 成長過程を原子レベルで明らかにした。AlN の MOVPE 成長においては原料ガスとして Al(CH3)3（TriMethyl Aluminum: TMAl）と NH3、キャリアガスとして N2 と H2 が用いられるが今回は MOVPE 成長中の TMAl の分解過程を量子論的に明らかにした。
今年度検討した TMAl 分解の経路を図 28に示す。図 28 の 29 個の反応過程に対して実際に第一原理計算を行って反応の自由エ
ネルギーを計算し、さらに NEB 法と String法による解析によって気相反応の活性化エネルギーと反応速度定数を求めた。計算で得られた TMAl の初期分解反応の活性化エ
ネルギーを図 29 に示す。図 29 からわかる
ように TMAl は初期は H2 ではなく NH3 と反応して Al(CH3)2(NH2)を生成する反応が

図 28：今回考慮した TMAl の気相中の分解反応支配的になることがわかる。詳しい計算の結果、GaN と同様の成長条件では TMAl は主に NH3 とのみ反応し、①Al(CH3)3 ＋ NH3 Al(CH3)2(NH2) ＋ CH4、②Al(CH3)2(NH2) + NH3  Al(CH3)(NH2)2 ＋ CH4、③Al(CH3)(NH2)2 + NH3
 Al(NH2)3 ＋ CH4 のように反応が進みアダクト系分子である Al(NH2)3 が最終生成分子として形成- 20 -

(a)(b)

図 29：TMAl の初期反応の模式図。(a)TMAl と NH3 の反応。(b)TMAl と H2 の反応。」

されることが明らかとなった。このように Al-N 結合を持つ Al(NH2)3 が AlN 表面に吸着することは、
表面の Al 液滴と気相中の NH3 が反応して成長すると一般に考えられているエピタキシャル成長の描像とは一致しない。実際に Eyring の式で決定した速度定数を用いて気相中の各化学種のモル分率の時間依存性をレート方程式

𝜕𝜕(𝜌𝜌𝑌𝑌𝛼𝛼 )𝜕𝜕𝜕𝜕= 𝜔𝜔𝛼𝛼 を解いて決定した。ここで、𝜌𝜌は密度、𝑌𝑌𝛼𝛼 は化学種𝛼𝛼の質量分率、

𝜔𝜔𝛼𝛼 は化学種𝛼𝛼の正味の生成速度である。その結果、GaN と同様の成長条件では、図 30(a)に示すようにアダクト系の分子である Al(NH2)3 が主な最終生成物として AlN 基板に飛来することがわかる。ところが、成長条件において NH3 原料ガスの割合(V/III 比)を GaN の 1/100 程度にし、キャリアガスである H2 ガスの割合を 4 倍の 0.8 にすると 0.1 秒程度のときには AlH3 が最終生成物となっているが、１秒程度経過すると AlNH2 や AlH が最終生成物になることがわかる(図 30(b))。実際の成長

(a)
TMAl

装置の形状と原料ガスの流入スピードを考慮
すると原料ガスは約 0.1 秒で AlN 基板に到達

Al(NH2)3するので（原料ガスの流入速度が約 1m/s、流入口から基板までの距離が 10cm 程度であAlH(NH2)2る）、この条件では基板近傍では AlH3 が最終生成物として AlN 基板に飛来することが予想される。このように、AlN の MOVPE 成長の場合には GaN に比べて V/III 比を 1/100 程度にして H2 キャリアガスのモル分率を増加さ

(b)
TMAl

せることで、Al-N 結合を含まない AlH3 が最終生成物として AlN 基板に飛来するようにな

AlH2NH2

り、成長表面に Al 液滴が形成されるというエ

AlH3

ピタキシャル成長の常識と一致する結晶成長

AlNH2
AlH

が可能になる。この結果は現実の AlN のMOVPE 成長で V/III 比を 1/100 にして H2 の分圧を大きくしていることと符合するも のであり、GaN に続いて AlN の MOVPE 成長の素過程を明らかにすることができた。

図 30：MOVPE 成長中の各化学種のモル分率の時間変化。(a)H2 のモル分率 0.2。(b)H2 のモル分率 0.8。

- 21 -

③ GaN MOVPE における気相反応のデータ同化 [九大計算グループ ＋ 名大計算／実験グループ]
名大計算グループによって GaN MOVPE の気相反応モデルが提案されており第一原理的に速度定数が計算されている。その結果示された主要反応経路は、名大実験グループで実施された重水素を供給した条件での質量分析（TOF-MS）実験の結果からもサポートされており、定性的な妥当性が示されている。今回、この反応モデルを定量的結晶成長シミュレータに応用するために、速度定数をチューニングするデータ同化を行った。
単に CH4 の TOF-MS データを正確に再現するように
チューニングを行うと、過適合することが容易に想像される。そこで、第一原理計算による活性化エネルギー

図 31：CH4 濃度の実験値と予測値の比較。

比を可能な限り維持するように多目的最適化を行った。
その結果、多数の速度定数セットの候補を得た。図 31 はそれらの予測性能を示しており、y=x の直線に近いプロットを与えた速度定数セットほど正確に予測している。チューニング前の速度定数セッ
トによる予測（図 31 中の×）から大きく改善していることが分かる。赤から紫までのマーカーの色は、
実験データをより重視するものから第一原理計算データをより重視するものまでのグラデーションに対応している。
こ れ ら の 速 度 定 数 を 用 い た 反 応 シ ミ ュ レ ー シ ョ ン を よ り 詳 細 に 見 て い く 。 図 32 は 左 図 がTMG+H2DMG-H+CH4、右図が TMG+NH3DMG-NH2+CH4 についての反応速度の時間変化を示している。実験データをより重視したもの（赤線）は、H2 による TMG 分解が支配的である結果を与えている。一方、第一原理計算データをより重視したもの（紫線）は、NH3 による TMG 分解が支配的である結果を与えており、これまでに支持されている主要反応経路を再現することができる。図 31 の紫マーカーを見ても十分に定量的な再現が可能と言える。以上のアプローチによって、定量的な予測性
能、第一原理的な基礎づけ、実験的・理論的に支持されている反応経路の 3 点を両立できる速度定数パラメータを得ることができた。

図 32：反応速度の時間変化。
（左図は TMG+H2DMG-H+CH4、右図は TMG+NH3DMGNH2+CH4）- 22 -

(4) GaN エピタキシャル成長素過程の解明：TOF-MS を用いた MOVPE 気相反応解析[名古屋大学実験グループ]
本テーマではこれまで GaN の MOVPE 成長における気相反応、すなわちトリメチルガリウム(trimethylgallium, TMGa, Ga(CH3)3)、トリエチルガリウム(triethylgallium, TEGa, Ga(C2H5)3)と N2、H2、アンモニア(NH3)との分解・反応について実験的な検証を行ってきた。今後は、AlN、InN など他の窒化物半導体やこれらの混晶に関する検討を進めていく計画である。AlN の Al 原料にはトリメチル
アルミニウム(trimethylaluminum, TMAl, Al(CH3)3)、InN の In 原料には TMIn が用いられる。特にAlN あるいは AlGaN の成長においては、NH3 との結合が強い TMAl を原料とするために、成長反応経路においてアダクト分子の形成、分解が主要な役割を果たすことが、理論計算および実験において示唆されてる[1-4]。窒化物半導体 MOVPE のあらゆる系において主要な N 源である NH3 の反応性について、正しく理解しておくことがより重要となる。また、前年度の検討においては、TMGa/H2/NH3 系における TMGa 分解・反応の解析を行うために、H2 の替わり重水素(deuterium, D2)をキャリアガスとして利用する同位体追跡法を用いる検討を行った。
この中で、NH3 がこの実験系において分解する 900℃よりも低い、400℃付近から D2 との交換反応を呈していることが明らかになった。これは、900℃程度まで安定であると考えられていた NH3 が、400℃において反応活性な状態となっていることを示唆しており、窒化物半導体 MOVPE の反応メカニズムを解析する上で、極めて重要な現象である。そこで今年
度は、D2 を用いた NH3 の反応活性化について、詳細な検討を行った。

図 33: 実験に用いた簡易 MOVPE 炉と分析システム。

- 23 -

図 34:

各種キャリアガス雰囲気での NH3 分解。

図 35: D2 キャリアガス中の NH3 反応。

本研究にあたっては前年度同様に図 33 に示す石英管状炉による簡易 MOVPE とマルチターン方式を採用した小型・高分解能な飛行時間型質量分析(time-of-flight mass spectrometry, TOF-MS)を組み合わせた実験系を用いた。主な実験条件は総流量 1000sccm、100kPa で一定とし、電気炉の設定温度を室温から 1200℃まで変化させている。石英反応管の内径は 26mm であり、室温における流速は0.02m/sec である。この流速は通常の MOVPE の成膜条件である 1m/sec に対して非常に遅い条件であり、本実験はガスのサンプリングは加熱領域端から 20cm の場所で行うことで、およそ平衡に近い状態で観察しているものと考える。実験は定めたガスフロー条件を維持し、ヒーター温度を段階的に挙げながら TOF-MS でマススペクトルを取得する。その後、着目する分子種のピーク強度の温度依存性をグラフ化することで解析を行っている。また、本実験系においては、NH3 自体のマススペクトルにおける信号強度は強すぎて、他の反応種と同時に観察するには適していないため、特に説明が無い場
合には主なフラグメント種である NH2 を NH3 と見立てて観察し、本文記載上はこれを NH3 と表現する。
まず、従来考えられている NH3 の分解反応を確認するため、MOVPE で用いられるキャリアガスである H2 と N2 に加え、N2 に変わる不活性ガス(inert gas)として He、Ar、をキャリアガスとして実験を行った。NH3 とキャリアガスの流量は それぞれ 50sccm と 950sccm である。図 34 にそれぞれのキャリア
ガス雰囲気で観察した NH3 強度の温度依存性を示す。変化量を比較しやすいように縦軸強度を室温の強度で規格化している。大別して見ると、NH3 強度は不活性ガスである N2、He、Ar 中では温度上昇と共に徐々に低下し、1150℃では 50%程度に達する。完全な平衡状態においては NH3 は 300℃程度で N2と H2 に分解するが[5]、その反応速度が非常に小さいため、低温域における緩やかな低下が見られていると考えられる。ガス種による違いは NH3 の熱拡散のしやすさ、気相熱伝導率の差異によるものと考えられる。特に He は N2、Ar に比べて衝突断面積が小さく、NH3 が流動しやすいため、同じ温度環境においても分解が促進されると考えられる。一方、H2 をキャリアガスとした場合には、他のガスに比べて高温域での分解が進みにくい事が解る。H2 分子は He とほぼ同じ衝突断面積であるにもかかわらず最も差異が大きい。これは H2 が反応活性であり、NH3 の N 原子から H 原子が脱離するのを抑制する効果があることは明白である。このように、NH3 の分解は温度と H 原子の存在によって影響を受ける。
H2 雰囲気における NH3 の分解率は 1150℃においても 20%程度であり、前述のように本実験系は非常に流速が遅い条件であるため、高流速条件においては気相の温度はより低い方にシフトし、NH3 の分解- 24 -

率はせいぜい数%であると考えられる。このことから、従来 NH3 は GaN の成長条件である 1000℃以上においても非常に安定であるとされている。
次に、キャリアガスとして H2 の代わりに同位体である D2 を用いた実験行った。NH3 と D2 の流量はそれぞれ 65sccm と 935sccm である。図 35 に NH3 および D2 との反応種強度の温度依存性を示す。縦軸は指数である。700℃程度までは NH3 の強度は安定しており、約 1/20 の強度で NH2D、約 1/200 の強度で
NHD2 が観察された。これは NH3 の H が D2 の D と置換して出来たものであるが、強度はほぼ一定あり、
温度依存性が見られない。700℃を越えると NH2D、NHD2 さらには ND3 の強度が急速に増大し、NH3 はそれに対応して低下している。NH2D は 950℃、NHD2 は 1030℃付近でピークとなり、減少に転じ、1150℃
においてはほとんどが ND3 となり、ほとんどの NH3 が周囲の D2 と置換してしまうことがわかる。この1150℃における ND3 の強度は室温における NH3 強度の半分程度であり、ここには示していないフラグ
メント成分である ND4 も含めると、総数としてほとんどのアンモニア種は分解せずに置換されたのみである。ここで重要なことは、図 34 において H2 雰囲気では緩やかに分解が進んでいるだけのように見える NH3 が、700℃以上では反応活性となり、H-D 交換反応が顕著となることである。これは、実際の結晶成長に用いられる高流速条件においても、NH3 が反応活性となっており、そのことが気相反応
あるいは表面反応においても影響を及ぼしているであろうと考えられる。

図 36：D2 キャリアガス中 TMGa 供給下での NH3 反応。

ここで改めて昨年の成果である TMGa/H2/NH3 および TMGa/D2/NH3 の実験結果を図 36 に示す。キャリ
アガスが H2 である場合は図 36(a)に示されるように、NH3 は 900℃当たりまで強度変化がなく、900℃
以上で分解が進んでいる。一方キャリアガスが D2 の場合はアンモニア種の総数には大きな変化がな
いが、H-D 交換反応は 400℃付近から発生しており、700℃にはほとんどが ND3 に置換されていることがわかる。キャリアガスを D2 にして観察すると、図 36(b)に示すように分解よりも低い温度で反応活性になり、交換反応が観察できることは図 35 の結果と同じではある。しかしながら、そこには 300℃
程度の大きな乖離が生じている。この活性化温度の違いは、TMGa の供給有無によるものであることは明らかである。従来 NH3 の分解においては、金属の触媒効果により活性化されることはよく知られている[5-7]。このため、TMGa の Ga が触媒として作用することは十分考え得ることである。しかしながら、気相中の TMGa 分子が直接触媒効果を呈するのか、実験中に上流で TMGa が分解し、石英管の炉壁に積層したものが触媒となっているかなど、不明な点が多い。そこで次のような切り分け実験を行っ- 25 -

た。

図 37: 触媒効果確認のための Ga 化合物設置外観。

実験は、反応炉の上流側に 3 種類、4 条件の Ga を含む化合物を設置し、その存在が NH3/D2 の反応・分解温度依存性にどう影響を与えるかを比較した。図 37 に各実験の外観写真を示す。実験①は図36(b)と同様に室温から TMGa を供給しながら分析を行うが、300℃で TMGa の供給を停止し、その後も昇温しながら測定を続ける。この場合、300℃までに上流で分解・反応した TMGa の一部が炉壁に残留していることになる。写真では解りにくいが、Ga と思われる白い堆積物が薄く付着している。実験②は上流に設置した石英ボートに金属 Ga を入れたものである。写真は室温の状態であるが、Ga は融点が 28℃であるため、実験中は融解して一つの液滴としてボート内に広がる。尚、Ga の蒸気圧は極めて低いため、1300℃以下では全く蒸発しないと考えてよい。実験③、実験④はサファイア基板上に
MOVPE で 2µm 成長した GaN ウエハを短冊状に切断したもので、表面は GaN (0001)面である。実験③は20mm×100mm 程度の切片を 1 枚、実験④は 4 枚を配置している。これらの条件で同様に D2/NH3 の分解・反応温度依存性実験を行った。図 38 に実験から得られた NH3 および ND3 強度の温度依存性を示す。実験①～④全ての条件において、NH3 強度低下、ND3 強度増大温度が何も設置しない場合に比べて低温化しており、Ga あるいはその化合部による触媒効果が確認された。更に NH3 強度低下温度で切り分けると、①TMGa300℃まで供給と②金属 Ga は同じ温度で 600℃、③④GaN ウエハでは共に 500℃となった。
また、図 36(b)の TMGa を高温まで連続供給した場合には、NH3 の強度低下は 400℃であった。このことから、同じ Ga を含む物質であっても、その状態によって効果が現れる温度が異なり、触媒機構において何らかの差異があるものと考えられる。①TMGa300℃まで供給時の NH3 強度低下温度は Ga 金属と同じであることから、TMGa 分解によって付着した白い堆積物は、推定通り金属 Ga に近い状態であると考えられる。実験③④での GaN ウエハを設置した場合、NH3 強度低下温度はほぼ同じであるが、設置枚数が多い④の方が分解が一気に進む一方で、③の場合には高温域で緩やかな低下がみられた。
この違いは触媒となる GaN の面積によるもので、③の場合は面積が十分でないために未活性な成分が残るためであると考えられる。一連の実験の中で、TMGa を高温まで連続供給した場合が最も NH3 強度低下温度、すなわち NH3 の反応活性化温度が低い結果となった。このことは、触媒効果が炉壁に付着した固相(または液相)の Ga 化合物だけではなく、気相中の TMGa 分子(およびその分解種)が直接的に- 26 -

触媒効果を発揮している可能性を示唆するものである。

図 38: 各種 Ga 化合物設置下での NH3 反応。

GaN の MOVPE 成長において、V 族元素すなわち N 原料である NH3 の反応性を明らかにするため、キャリアガス H2 の代わりに同位体の重水素(D2)を用いて、同位体追跡による質量分析をおこなった。従来の NH3 の反応性は NH3 が N2 と H2 に分解する、分解反応により議論していたが、分解する温度よりも遥かに低い温度で、NH3 分子の H と D2 分子の D が交換反応することを明らかにし、NH3 が従来よりも低い温度で反応活性な状態になっていることを見出した。また、TMGa 供給、金属 Ga、GaN 結晶によってNH3 の活性化温度を低温化する触媒効果を確認した。この NH3 活性化と触媒効果は、研究項目(3)の成長素過程解明量子論計算において、気相反応や GaN 結晶表面における NH3 と Ga 原子の相互作用の理解にも貢献できる成果である。
参考文献[1] T.G. Mihopoulos et al., J. Cryst. Growth 195, 733Ð739 (1998).
[2] J. Stellmach et al., J. Cryst. Growth 315, 229 (2011).
[3] K. Matsumoto et al., J. Cryst. Growth 272, 360 (2004).
[4] G.S. Huang et al., Mater. Sci. Eng. B 136, 29 (2007).
[5] S.F. Yin et al., Applied Catalysis A: General 277, 1 (2004).
[6] S.S. Liu et al., J. Electrochem. Soc. 125, 1161 (1978).
[7] K. Sasamono et al., Phys. Status Solidi C 8, 2092 (2011).
[8] K. Sasamono et al., J. Cryst. Growth 314, 62 (2011).
(5) プロジェクトの総合的推進
プロジェクト推進のための打ち合わせはオンラインおよび対面で適宜開催した。とくに名古屋大学実
験グループと計算グループの間では、日常的に議論が行われた。また本事業参加のメンバー全員が集結
して、2021 年 7 月 15 日にオンラインで全体会議を行った。詳細と発表資料は https://fugaku-- 27 -

semicon.jp/achievements/20210715.html に掲載中である。
また「富岳」成果創出加速プログラムの他の課題との連携をはかるため、材料系 7 課題の合同研究会（参加人数 139 名）を 2021 年 12 月 8 日-9 日にオンラインで開催した。シミュレーション技術に関する有益な情報が交換された。詳細と発表資料は https://fugaku-semicon.jp/jointmtg21/ に掲載されている。
また、本プロジェクト主催の国際会議として、2022 年 1 月 11 日-13 日に名古屋国際会議場およびオンラインで、The 2nd International Symposium on Wide Gap Semiconductor Growth, Process andDeviceSimulation(ISWSGPDs) が 開 催 さ れ た 。 詳 細 は 、 https://sites.google.

com/view/iwsingularity3iswgpds2/home に掲載中である。
また、本プロジェクトのウェブページも日英 2 言語ですでに作成済みであり、多くの訪問者を数えている（https://fugaku-semicon.jp/）。
(6)新たな密度汎関数法（DFT）計算手法の開発
最後に、2-1 の実施項目には記載がないが、上記研究項目(1), (2), (3)に共通な計算技術基盤として、新たに、電子軌道の計算を必要としないオーダーN・DFT 法を、ニューラルネットワークを活用して開発し、その有効性を実証したので、以下に報告する。また、将来の量子コンピューティングの発展を見据えた新たなアルゴリズムの開発も報告する。
①

ニューラルネットワークを活用した新たな DFT 計算手法の開発[名古屋大学計算グループ]
DFT の基本定理は、相互作用しあっている電子系（つまり実際の物質）の全エネルギー E は電子

密度 ρ (r ) の汎関数として正確に記述できる（Hohenberg-Kohn の定理）というものである。したがって、その汎関数の形がわかれば E[ ρ ] の ρ (r ) に対するオイラー方程式を解くことにより、 ρ (r ) が求まり、全エネルギーが計算でき、ターゲット物質の物性が議論できる。しかし E[ ρ ] の正確な形はわからない（それは量子多体問題を解くのと同義である）。そこで全エネルギーを独立な仮想電子系の運動エネルギーとその残りの部分に分けるスキームが始められた（Kohn-Sham スキーム）。このス
キームでは仮想電子系の運動エネルギーは各電子の軌道 ϕi (r )（Kohn-Sham 軌道： i = 1...N ）であら
わされ、もともとの ρ に対するオイラー方程式は ϕi に対するターゲットサイズ N に比例する数の方程式群（Kohn-Sham 方程式）に変換される。この KS スキームは、多くの物質群に対して大きな成功をおさめ、計算物質科学の隆盛を支えている。しかし N 個の軌道を導入したおかげで、Kohn-Sham
方程式は非線形に絡み合った方程式群となり、その解法の計算コストはオーダー N 3 となっている。
しかし DFT の基本に立ち返れば、運動エネルギーの、 ρ についての汎関数の形がわかれば、解くべきオイラー方程式は一つであり、計算コストはオーダー N となり、飛躍的な計算効率の上昇が見込まれる。この軌道を用いないスキームを Orbital-Free DFT （OFDFT）とよぶ。今回、この運動エネル
ギー汎関数の形をニューラルネットワーク(NN)を用いて明らかにし、新たな NN-OFDFT 法を開発した。
運動エネルギー汎関数 T [ ρ ] の形としては以下のようなものを考えた。

T [ρ ] = ∫τTF(r ) F ( s 2 , q; W )dr- 28 -

ここで、τ TF は一様電子密度の系に対する運動エネル
ギー密度、F は実際の系での、一様電子密度の場合の運動エネルギーからのずれを表す enhancement 因子である。 F の形としては様々なものを考えることが可能だが、ここでは電子密度とその 2 次微分までの関数で表されると仮定する（2 次の密度勾配近似）。す2 1/3

な わ ち s = ∇ρ [2(3π )ρ 4/3 ] 、 お よ び

q = ∇ 2 ρ [4(3 π 2 )3/2 ρ 5/3 ] である。W が今回導入したNN の重みパラメータ W jk(l ) の集合であり、 l − 1 層の

図 39: ニューラルネットワークの模式図。
hidden layer の数が 2 の場合。q と s が記述子であり、enhancement 因子 F が出力である。実際の計算では 3 個の hidden layer、
layer あたり 5 個のニューロンを用いた。

k 番目のニューロンと l 層の j 番目のニューロンとを結びつけるシナプスと言える。
（図 39）この重みパラメータを training により教育し、最適な運動エネル
ギー汎関数を求める。
Training では以下のコスト関数 L を最小化する。

1 Nt 1  δ T (rn ) δ T KS (rn ) 
=
−
L

∑ 
δρ 
Nt n=1 2  δρ2

すなわち運動エネルギー汎関数の 1 次の微分（1 次汎関数微分）の値が実空間の各メッシュ点 rn で、 KohnSham スキームで得られる値 δ T KS δρ と可能な限り一{ }致するように W jk(l ) を定める。実際の最適化は natural
gradient descent 法で行った。Training はダイヤモン
ド構造の炭素の計算データを用いて行った。そこで得られた運動エネルギー汎関数が他の物質に対して有効
かどうかが問題である。
図 40 はシリコン結晶に対して、今回の汎関数で得られた電子密度分布(NN)と Kohn-Sham スキームで得ら

図 40：ダイヤモンド構造 Si の(111)方向に沿った電子密度。ニューラルネットワークで得られた汎関数（NN）と Kohn-Sham スキームで得られた結果（KSDFT）、および過去の汎関数による結果。(b)は KSDFT とのずれがプロットされている。

れた密度分布(KSDFT)を比較したものである。図中の TF(1/5)vW、LKT、PGSL0.25 等は過去の近似理論で得られた汎関数を用いた OFDFT 計算結果である。今回の NN 汎関数が最も Kohn-Sham スキームの結果に近いことが見て取れる。Si 以外にも金属からイオン結晶に至る 24 種の異なる物質群、結晶多型（diamond, graphene, fcc-Si, beta-tin-Si, SiC, bcc-Li, fcc-Al, fcc-Cu, bcc-Na, NaCl および 7 種の原子、6 種の分子）に対しても電子密度が計算され、その各メッシュ点での値の Kohn-Sham スキー
ムからのずれの RMSE（root mean square error）は、従来の最良の汎関数のずれよりも、39%減少することが判明した。

- 29 -

表 2：4 種の金属、イオン結晶、グラフェンの構造定数。MARE は他の 20 種の物質に対する構造定数と合わせて Mean Absolute Relative Error を示したもの。

また構造定数（固体の場合の格子定数、体積弾性率、分子の原子間距離等）も計算され、KohnSham スキームで得られた値と良好な一致を示している。一例として、表 2 に fcc-Cu、 bcc-Na、
NaCl、graphene の例を示す。また他のすべての物質群に対する構造定数計算結果と合わせて、すべての物質群に対しての Kohn-Sham スキームによる結果からのずれの MARE (Mean Absolute
Relative Error)も表 2 に示してある。過去の汎関数に比べての優位性は顕著である。ダイヤモンドのデータで学ばせたニューラルネットワークが、このように性質の異なる広範な物質群に対して良好な結果を与えることは、ひとつの驚きである。
図 41 は、今回 開発した OFDFT コード と

図 41: OFDFT（青）と KSDFT（赤）との計算時間の比較。テスト系は SiC の超格子で 576 原子から 4704 原子までの結果。(b)が SCF 全体計算の結果であり、(a)では各計算部分に分割しての計測結果。t1: 部分空間対角化、t2: 共役勾
配法計算とベクトル直交化、t3: 電子密度の更新。

KSDFT で最速とみなされる RSDFT コードの計算時間比較である。対象は SiC 結晶を選んだ。今回の NN-OFDFT コードはきれいにオーダーN の計算時間スケーリングを示している。ミクロンスケー
ルの超大規模計算への応用が期待できる。
以上の成果は、Physical Review Research 誌に公表された。
②

量子コンピューティングでのハイブリッドアルゴリズムの開発[東京工業大学計算グループ + Quemix 計算グループ]
材料科学を飛躍的に発展させ得る、量子コンピュータのアルゴリズムの開発にも着手した。材料科学において昔から知られた虚時間発展法に注目し、量子コンピュータ上でのアルゴリズムを開発した。特に、古典コンピュータとの良いとこ取りをする、ハイブリッド方式量子アルゴリズムの開発を行った。虚時間発展法を量子コンピュータ上で実装する際に、一番の課題は非ユニタリ演算で- 30 -

ある虚時間発展演算子をいかにしてユニタリゲート操作のみの量子コンピュータ上で実装するか、である。本年度は、虚時間発展演算子が作用されたベクトルをよく近似するようにユニタリ演算子
を古典コンピュータにより求めて、量子コンピュータ上でゲート操作として実行を行うハイブリッド方式量子アルゴリズムの開発に成功した。特に、ユニタリゲート操作として近似する際に部分ハ
ミルトニアンの持つ局所性を拡張し、非局所化させることにより、量子回路深さを大幅に低減させると同時に、近似の精度を改善できることを明らかにした。図 42 には、開発した非局所近似(NLA)の適用結果を示す。図 42 は、非局所近似虚時間発展法を 10 サイトの 3 正則マックスカット問題に適用した結果である。この問題において、厳密解は-12 であり、虚時間発展ステップ 10 において厳密解に対して 2 次の非局所近似（NLA-D2）にて 95%の精度、3 次の非局所近似（NLA-D3）にて100%の精度を出していることがわかる。一方、古典アルゴリズムを用いた際の近似精度では 93%であり、量子虚時間アルゴリズムを用いた計算では古典アルゴリズムの精度限界を超え得る可能性を見出すことができた。また、局所近似の制限範囲内では、2 次の局所近似（LA-D2）にては非単調減
少な振る舞いが見られ、非局所近似を用いることで精度が改善していることがわかる。得られた成果は、npj quantum information 7, 85(2021)にて報告を行なった。

図 42: 虚時間発展法を 10 サイトの 3 正則グラフへと適用した際の量子コンピュータのシミュレータ結果。横軸は虚時間ステップ、縦軸はエネ
ルギーを示す。局所近似（LA）、非局所近似(NLA)、拡張局所近似(eLA)の比較を行う。

２－３．活動（研究会の活動等）
以下のシンポジウム、ワークショップを主催した。
1.

The 2nd International Symposium on Wide Gap Semiconductor Growth, Process and Device Simulation(ISWSGPDs), Nagoya Congress Center (Hybrid format), 2022.01.11-13.

2.

第 39 回コンピューテーショナルマテリアルズデザインワークショップ, アドバンストコース, 2021
年 9 月 6 日～2021 年 9 月 10 日, オンライン （大阪大学）

3.

第 40 回コンピューテーショナルマテリアルズデザインワークショップ, スパコンコース, 2022 年 2月 21 日～2022 年 2 月 25 日, オンライン. （大阪大学）- 31 -

２－４．実施体制
業務項目

担当機関

（１）量子論第一原理計算に〒464-8601

よる表面・界面原子機構解明

愛知県名古屋市千種区不老町と実験との比較国立大学法人東海国立大学機構

担当責任者

押山淳

名古屋大学
〒152-8550
東京都目黒区大岡山 2-12-1
国立大学法人東京工業大学

松下雄一郎〒657-8501
兵庫県神戸市灘区六甲台町 1-1
国立大学法人神戸大学（２） 統合デバイスシミュ〒565-0871

レータによる省エネルギー

大阪府吹田市山田丘 2-1

MOSFET の特性解明と実デバ

国立大学法人大阪大学

イスとの比較

〒464-8601

小野倫也

森伸也

愛知県名古屋市千種区不老町
国立大学法人東海国立大学機構
名古屋大学

押山淳

（３） エピタキシャル成長〒464-8601

統合シミュレーション技術の愛知県名古屋市千種区不老町

開発と有用性の検証

国立大学法人東海国立大学機構
名古屋大学

白石賢二

〒816-8580
福岡県春日市春日公園 6-1
国立大学法人九州大学（４）GaN エピタキシャル成〒464-8601

長素過程の解明：TOF-MS を愛知県名古屋市千種区不老町

用いた MOVPE 気相反応解析

国立大学法人東海国立大学機構
名古屋大学

寒川義裕

天野浩

（５）プロジェクトの総合的〒464-8601

推進愛知県名古屋市千種区不老町
国立大学法人東海国立大学機構
名古屋大学

押山淳

- 32 -

- 33 -

別添１ 学会等発表実績
発表論文（査読あり）
1.

K. Seino and A. Oshiyama, “Microscopic mechanism of adatom diffusion on stepped SiC
surfaces revealed by first-principles calculations" Appl. Surf. Sci. 561, 149927 (2021).

2.

K. Chokawa, K. Shiraishi and A. Oshiyama, “Defect-free interface between amorphous(Al2O3)1-x(SiO2)x revealed by first-principle molecular dynamics simulations" Appl. Phys.
Lett. 119, 011602 (2021)

3.

F. Imoto, M. Imada and A. Oshiyama, “Order-N orbital-free density-functional calculations
with machine learning of functional derivatives for semiconductors and metals" Phys. Rev.
Research 3, 033198 (2021).

4.

A. Kusaba, Y. Kangawa, T. Kuboyama and A. Oshiyama, “Exploration of a large-scalereconstructed structure on GaN(0001) surface by Bayesian optimization" Appl. Phys. Lett. 120,021602 (2022).

5.

A. Kusaba, R. M. Kurniawan, P. Kempisty, Y. Kangawa, “DFT modeling of unintentional oxygen
incorporation enhanced by magnesium in GaN(0001) and AlN(0001) growth surfaces during
metal-organic vapor-phase epitaxy”, Phys. Status Solidi B, 2022, 2100430. (in Press)

6.

N. Mori, F. Hashimoto, T. Mishima, and H. Tanaka, "Analytical models for inter-layer tunnelingin two-dimensional materials" Jpn. J. Appl. Phys. 61, SC1022 (2022).

7.

J. Okada, F. Hashimoto, and N. Mori, "Equivalent model for band-to-band tunneling simulationof direct-gap III–V semiconductor nanowires" Jpn. J. Appl. Phys. 60, 091002 (2021).

8.

N. Mori, A. Komada, and G. Mil'nikov, "Comparison of linear and quadratic dispersion models
for phonon transport in one-dimensional mass-disordered systems" APL Materials, 9, 08112(2021).

9.

R. Negishi, K. Yamamoto, H. Tanaka, S. A. Mojtahedzadeh, N. Mori, Y. Kobayashi, "Crossover
point of the field effect transistor and interconnect applications in turbostratic multilayer
graphene nanoribbon channel" Sci. Rep. 11, 10206 (2021).

10. M. Uemoto, N. Komatsu, Y. Egami and T. Ono, "First-Principles Study on Structure andAnisotropy of High N-atom Density Layer in 4H-SiC" J. Phys. Soc. Jpn. 90, 124713 (2021).
11. M. Tsunasaki, T. Ono, and M. Uemoto, "Theoretical investigation of vacancy related defects at
4H-SiC(000-1)/SiO2 interface after wet oxidation" Jpn. J. Appl. Phys. 61, SH1001 (2022).
招待講演
1.

A. Oshiyama, “Computics Approach to Dislocation-Impurity Complexes and Interface
Characteristics of GaN Devices" The 31st Int Conf Defects in Semiconducotors (Oslo, Norway,
July 26-30, 2021, online).

2.

(Plenary) A. Oshiyama, “Computics Approach in Science of Power Electronics" 1st Int Conf on
Computational Science and Data Analytics" COMDATA 2021 (November 21-24, 2021, Kuala
Lumpur, Malaysia).

- 34 -

3.

(Keynote) Akira Kusaba, Application of Machine Learning Methods to More Quantitative GaN
MOVPE Modeling, The 2nd International Symposium on Wide Gap Semiconductor Growth,
Process and Device Simulation, Nagoya (hybrid format), 2022.01.11-13.

4.

T. Ono, "DFT study on defects at SiC MOS interface: electronic structure and formation
mechanism" The 3rd International Workshop on Materials Science and Advanced Electronics
Created by Singularity (IWSingularity 2022)/The 2nd International Symposium on Wide Gap
Semiconductor Growth, Process and Device Simulation (ISWGPDs 2022) (January 11-13, 2022,
Nagoya, Japan).

5.

Kenji Shiraishi, First Principles Studies of Impurity Dislocation Complexes in GaN, The 11th
International Conference on Processing & Manufacturing of Advanced Materials, 2021.06.1-5.
Online.

6. Kenji Shiraishi, First Principles Studies on the Atomistic Processes of GaN Metal Organic
Vapor Phase Epitaxy (MOVPE), 5th Edition of Applied Science, Engineering and Technology
Virtual, 2021.12.13-14, online.
7. T. Akiyama, T. Shimizu, T. Ito, H. Kageshima, and K. Shiraishi, Ab initio-based approach forreaction process at 4H-SiC/SiO2 interfaces, The 2nd International Symposium on Wide Gap
Semiconductor Growth, Process and Device Simulation, Nagoya (hybrid format), 2022.01.1113.
8. K. M. Bui, K. Shiraishi, and A. Oshiyama, Density Functional Theory study of Step Flow
Growth of Gallium Nitride, The 2nd International Symposium on Wide Gap Semiconductor
Growth, Process and Device Simulation, Nagoya (hybrid format), 2022.01.11-13.
9.

G. Mil'nikov and N. Mori, "Application of the R-matrix approach to large scale NEGF
simulations" Int Symp on Wide Gap Semiconductor Growth, Process and Device Simulation
2022 (January 11-12, 2022, online).

10. N. Mori, S. Nagamizo, H. Tanaka, and G. Mil'nikov, "Semiconductor device modeling andsimulation from an atomistic view" Int Sympo on Wide Gap Semiconductor Growth, Process
and Device Simulation 2022 (January 11-12, 2022, online).
11. Shugo Nitta, Ye Zheng, Daisuke Yahara, Yoshio Honda, and Hiroshi Amano, "High-resolution
vapor phase reaction analysis in nitride semiconductors MOVPE"

The 3rd International

Workshop on Materials Science and Advanced Electronics Created by Singularity(IWSingularity 2022)/The 2nd International Symposium on Wide Gap Semiconductor Growth,
Process and Device Simulation (ISWGPDs 2022) (January 11-13, 2022, Nagoya, Japan).
12. H. Tanaka and N. Mori, "Monte Carlo simulation of electron mobility in SiC MOSFETs" 2021
Int Meeting for Future of Electron Devices, Kansai (November 19, 2021, online).
13. （基調講演）白石賢二、本当のパワーデバイス材料をデータサイエンスで作るには？、複相機能
開拓拠点ワークショップ、東京（ハイブリッド）、2021.11.2.
14. 小川湧太郎，洗平昌晃，遠藤哲郎，白石賢二、Fe/MgO 界面への窒素不純物が磁気異方性と TMR に与える影響について、第 82 回応用物理学会秋季学術講演会、オンライン、2021.09.11-14.

- 35 -

15. 寒川義裕、窒化物半導体エピタキシーにおける不純物混入機構の理論解析、第 35 期コンピュータによる材料開発・物質設計を考える会（Computer aided Materials and Molecular Design
Forum: CAMM）フォーラム本例会、オンライン、2022.02.04.
16. 寒川義裕、窒化物半導体の薄膜成長におけるプロセスインフォマティクス 、結晶工学×インフォマティクスセミナー（第 26 回結晶工学セミナー、第 3 回インフォマティクス応用研究グルー
プ研究会）、オンライン、2021.11.16.
17. 寒川義裕、プロセス・インフォマティクスの新機軸：窒化物半導体 MOVPE における取り組み、日本学術振興会 R032 委員会第 2 回研究会、オンライン、2021.8.6.
18. 草場彰、久保山哲二、寒川義裕、結晶成長デジタルツイン―AI と計算科学からのアプローチ、
第 50 回結晶成長国内会議、オンライン、2021.10.27.
19. 田中 一，森 伸也，"界面準位を考慮した SiC MOS 反転層における Hall 移動度の理論解析"，
応用物理学会春季学術講演会，2022 年 3 月 23 日，オンライン.
20. 森伸也，星野 知輝，田中 一，美里劫 夏南，"パワーデバイスのモデリング・シミュレーショ

ン"，第 40 回電子材料シンポジウム，2021 年 10 月 12 日，オンライン.
国際会議発表1.

A. Kusaba, Y. Kangawa, Z. Ye, S. Nitta, K. Shiraishi, Tuning of Reaction Rate Constants for
Trimethylgallium Decomposition by Multiobjective Genetic Algorithm with High-Resolution
Mass Spectrometry Data, International Conference on Materials and Systems for
Sustainability 2021 (ICMaSS 2021), online, 2021.11.04.

2.

Y. Kangawa, A. Kusaba, P. Kempisty, Theoretical approach to unintentional oxygen doping
during MOVPE of GaN:Mg and AlN:Mg, International Conference on Materials and Systems
for Sustainability 2021 (ICMaSS 2021), online, 2021.11.04.

3.

A. Kusaba, P. Kempisty, Y. Kangawa, First-principles study of Mg and O co-doping
mechanism in the growth surface during GaN(0001) and AlN(0001) metalorganic vapor phase
epitaxy, 32nd IUPAP Conference on Computational Physics (CCP 2021), online, 2021.08.04.

4.

M. Tsunasaki, M. Uemoto, and T. Ono, "DFT Study on Defect Structures at SiC(000-1)/SiO2
Interface after Wet Oxidation" 2021 International Workshop on DIELECTRIC THIN FILMSFOR FUTURE ELECTRON DEVICES (November 14-16, 2021, Online).

5.

N. Komatsu, Y. Egami, T. Ono, and M. Uemoto "DFT study on structure and anisotropy ofhigh N-atom density layer in 4H-SiC" 2021 International Conference on Solid State Devices
and Materials (SSDM2021) (September 6-9, 2021, Sapporo, Japan).

6.

S. Nagamizo, H. Tanaka, and N. Mori, "Analysis of electronic states at SiC MOS interface
based on empirical pseudopotential method" 2021 Int Meeting for Future of Electron Devices,
Kansai (November 19, 2021, online).

7.

A. Komada and N. Mori, "Impact of anharmonic phonon-phonon scattering on phonontransport in one-dimensional system" 2021 Int Meeting for Future of Electron Devices,
Kansai (November 18, 2021, online).

- 36 -

8.

Y. Kunimoto and N. Mori, "Wigner Monte Carlo simulation of quantum superposition states"
2021 Int Meeting for Future of Electron Devices, Kansai (November 18, 2021, online).

9.

N. Mori, F. Hashimoto, T. Mishima, and H. Tanaka, "Analytical models for inter-layer tunnelingin two-dimensional materials" 2021 International Conference on Solid State Devices andMaterials (September 6-9, 2021, online).

10. S. Nagamizo, H.Tanaka, and N. Mori, "Electronic states in 4H-SiC MOS inversion layers
considering crystal structure using empirical pseudopotential method" Int Workshop on
Computational Nanotechnology (May 26-June 6, 2021, online).
11. J. Okada, F. Hashimoto, and N. Mori, "Equivalent model for tunneling simulation of direct-gap
semiconductor nanowires" Int Workshop on Computational Nanotechnology (May 26-June 6,
2021, online).
12. T. Kimura, K. Chokawa, A. Oshiyama and K. Shiraishi, Theoretical Study of H Adsorption on
Stepped SiC Surface during CVD Growth, 2021 International Conference on Solid State
Devices and Materials (SSDM2021), 2021.09.06-09, online.
13. N. Inoue, Y. Harashima, K. Chokawa, K. Shiraishi, and A. Oshiyama, First-principles
Calculations between Screw Dislocation with Mg, H Impurities on GaN, International
Conference on Materials and Systems for Sustainability 2021 (ICMaSS 2021), online,
2021.11.04-06.
14. F. Nanataki, M. Araidai, H. Kageshima, and K. Shiraishi, First-principles studies on theeffects of O atoms in the substrate on the oxidation of a vertical Si nanopillar, International
Conference on Materials and Systems for Sustainability 2021 (ICMaSS 2021), online,
2021.11.04-06.
15. R. Sakurai, K. Chokawa, A. Oshiyama, and K. Shiraishi, First Principles Studies on Atomic
and Electronic Structures of VGa-VN divacancies, International Conference on Materials andSystems for Sustainability 2021 (ICMaSS 2021), online, 2021.11.04-06.
16. Y. Ohata, M. Araidai, Y. Shibata, G. Hashiguchi, and K. Shiraishi, Effects of carbon atoms on
the reliability of the potassium-ion electret used in vibration-powered generators,
International Workshop on DIELECTRIC THIN FILMS FOR FUTURE ELECTRONDEVICES – SCIENCE AND TECHNOLOGY – (2021 IWDTF), online, 2021.11.14-16.
17. T. Kimura, K. Chokawa, A. Oshiyama, and K. Shiraishi, First-Principles Investigation of
Hydrogen Adsorption on Stepped SiC Surface during CVD Growth, The 2nd International
Symposium on Wide Gap Semiconductor Growth, Process and Device Simulation, Nagoya(hybrid format), 2022.01.11-13.
国内会議での一般講演・ポスター発表は省略します受賞- 37 -

1. 草場彰、第 19 回奨励賞、非平衡量子熱力学による GaN 気相成長プロセスの解明、日本結晶成長学会、2021.10.28.
2. 木村友哉，安田賞、CVD 成長環境下における SiC 微斜面への水素被覆の理論研究、第 27 回電子
デバイス界面テクノロジー研究会、2022.1.29.
3. 大畑慶記，Young Award, International Workshop on DIELECTRIC THIN FILMS FOR
FUTURE ELECTRON DEVICES – SCIENCE AND TECHNOLOGY – (2021 IWDTF),
2021.11.16.

- 38 -