修 士 論 文 の 和 文 要 旨 研究科・専攻 大学院 情報理工 学研究科 情報・ネットワーク工学専攻 博士前期課程 氏 名 平野 安彦 学籍番号 1631123 論 文 題 目 シミュレーションによるジャロシンスキー・守谷相互作用測定法の検討 要 旨 コンピュータ等の情報機器の記憶装置として現在主にDRAMが使われている。DRAMは情報の保持の ために電力の供給が必要であるという揮発性を有しており、情報機器の消費電力低減化を妨げる 一因となっている。電力の供給がなくても情報を失わない不揮発性メモリの一種としてMRAMがあ る。MRAMの実現により情報機器の消費電力の低減が期待されているが、実用化には高記憶密度化、 高速書き込み等の課題がある。これらの課題に対して「スカーミオン構造」と呼ばれる微細構造 や、「高速磁壁移動」等の、特異な磁性現象が生じる「ジャロシンスキー・守谷相互作用(DMI)」 が現在注目されている。これらの磁性現象により記憶セルの小型化や情報書き込み速度の高速化 についての報告が行われており、DMI及びDMIによる磁性現象についての研究が近年盛んに行われ ている。また、DMIの強度を示すDMI定数を測定する手法について様々な手法が提案されている。 DMI定数の測定には、直接磁性体の磁化の捻れを測定する方法や、DMIによる磁性現象を利用した 測定法等がある。本研究ではDMI定数の測定法について、新たに磁化反転実験による測定法を提案 し、マイクロマグネティックシミュレーションにて検討した。具体的には、反転磁界から求める 手法と、反転確率から求める手法の2つの手法である。反転磁界によるDMI定数測定法は、直径の 異なる複数の円盤状の磁性体に対して面直反転磁界及び面内反転磁界を求め、各磁界の差から測 定する手法である。提案手法を確認したところ、円盤直径の増加と共に面直及び面内反転磁界が 減少するが、DMI定数では面内反転磁界のみが変化することがわかった。面直反転磁界と面内反転 磁界の差を取ったところ、両者の差が0となる円盤直径がDMI定数により変化することがわかった。 これより反転磁界よりDMI定数を測定できることがわかった。反転確率によるDMI定数測定法は、 円盤状の磁性体に対して、磁界及び磁界パルス幅による反転確率の変化から測定する手法である。 シミュレーションで提案手法を確認したところ、パルス幅1 nsではDMI定数により反転確率が大き く変化するが、パルス幅が3 ns以上の長パルスではDMI定数による反転確率の変化が小さいことが わかった。また異なる磁界パルス幅による反転確率の変化量を「反転確率面積」とした場合、DMI 定数の増加と共に反転確率面積が単調に減少し、D=1.0 erg/cm2で約30%減少することがわかった。 これより反転確率よりDMI定数を測定できることがわかった。
電気通信大学情報理工学研究科 情報・ネットワーク工学専攻情報数理工学プログラム
平成
29
年度 修士論文
シミュレーションによる
ジャロシンスキー・守谷相互作用測定法の検討
平成 30 年 3 月 5 日提出情報数理工学プログラム
学籍番号
1631123
平野 安彦
指導教員 仲谷栄伸目 次
第 1 章 はじめに 5 1.1 研究の背景 . . . 5 1.2 目的 . . . 5 1.3 論文構成 . . . 6 第 2 章 基本事項 7 2.1 マイクロマグネティックシミュレーション . . . 7 2.2 磁気モーメント . . . 7 2.3 LLG方程式 . . . 8 2.4 ジャロシンスキー・守谷相互作用 (DMI) . . . 9 2.4.1 界面ジャロシンスキー・守谷相互作用 (界面 DMI) . . . 9 2.4.2 DMIの測定 . . . 9 2.5 実効磁界 . . . 9 2.5.1 異方性磁界 . . . 10 2.5.2 交換磁界 . . . 10 2.5.3 静磁界 . . . 10 2.5.4 DMI磁界 . . . 10 2.6 熱揺らぎ . . . 10 第 3 章 マイクロマグネティックシミュレーション 11 3.1 LLG方程式の数値解法 . . . 11 3.2 オイラー法 . . . 12 3.3 4次のルンゲ・クッタ法 . . . 13 3.4 計算領域の離散化 . . . 13 3.5 実効磁界の計算 . . . 14 3.5.1 異方性磁界 . . . 14 3.5.2 交換磁界 . . . 15 3.5.3 DMI磁界 . . . 16 3.5.4 境界条件 . . . 18 3.5.5 静磁界 . . . 19 3.6 静磁界の高速計算 (離散高速フーリエ変換) . . . 22 3.6.1 離散高速フーリエ変換を用いた Convolution 演算の高速計算 . . . 22 3.6.2 非周期構造での Convolution 演算 . . . 23 3.7 熱揺らぎ磁界における乱数生成 . . . 24第 4 章 シミュレーションモデル 26 4.1 マクロスピンモデル . . . 26 4.2 マイクロマグネティックモデル . . . 27 第 5 章 提案する DMI 定数の測定法 28 5.1 反転磁界による測定 . . . 28 5.2 反転確率による測定 . . . 28 第 6 章 反転磁界による DMI の測定法 30 6.1 面直反転磁界の計算 . . . 30 6.1.1 計算条件 . . . 30 6.1.2 計算結果 . . . 30 6.2 面内反転磁界の計算 . . . 31 6.2.1 計算条件 . . . 32 6.2.2 計算結果 . . . 32 6.3 面直、面内反転磁界の考察 . . . 34 6.3.1 実効異方性磁界との比較 . . . 34 6.3.2 DMI磁界による反転磁界への影響 . . . 34 6.3.3 外部磁界による面内磁化の変化 . . . 36 6.3.4 面直、面内反転磁界の差 . . . 36 6.4 まとめ . . . 38 第 7 章 反転確率による DMI の測定法 39 7.1 シミュレータの妥当性の調査 . . . 39 7.1.1 計算条件 . . . 39 7.1.2 計算結果 . . . 40 7.2 計算の安定性の調査 . . . 40 7.2.1 計算条件 . . . 40 7.2.2 計算結果 . . . 43 7.3 反転確率の計算 . . . 43 7.3.1 計算条件 . . . 43 7.3.2 計算結果 . . . 44 7.4 磁界パルス幅による反転確率変化の定量的比較 . . . 44 7.4.1 比較方法 . . . 45 7.4.2 計算結果 . . . 46 7.5 まとめ . . . 46 第 8 章 反転確率の外部磁界印加角度依存性 48 8.1 外部磁界印加角度による反転確率のシミュレーション . . . 48 8.1.1 計算条件 . . . 49 8.1.2 計算結果 . . . 49 8.2 外部磁界印加角度による DMI 測定法への影響 . . . 49 8.2.1 計算条件 . . . 52 8.2.2 計算結果 . . . 52
8.3 低確率での反転確率の傾き . . . 53 8.3.1 計算結果 . . . 53 8.4 まとめ . . . 56 第 9 章 まとめ 57 参考文献 60 発表リスト 61 謝辞 62
第
1
章 はじめに
この章では、本研究の背景、目的、本論文の構成について述べる。まず研究の背景では、記 憶装置の性能向上に関わる磁性現象の背景にあるジャロシンスキー・守谷相互作用及び DMI 定 数の測定法について述べる。目的では、本研究における目的について説明する。最後の論文構 成では、各章における内容について述べる。1.1
研究の背景
現在、コンピュータなどの情報機器の主記憶装置には、主に DRAM(Dynamic Random Access
Memory)が用いられている。DRAM は情報の記録を電荷により行うため、常に電源を供給し
なければ情報を記録できなくなる。この性質を揮発性という。DRAM のように揮発性を有する メモリは揮発性メモリと呼ばれる。これに対して、情報の記録に電源の供給が必要ない不揮発 性のメモリの一つに MRAM(Magnetoresistive Random Access Memory) がある。MRAM は、 情報の「0」「1」を記憶セルに用いる磁性体の磁化の向きで記録するため、記録の保持には電源 供給の必要が無い。そのため、DRAM に比べて消費電力の低減が見込まれる。しかし、現在製 品化されている MRAM は DRAM に比べて記憶容量が小さく、記憶容量増大のために記憶セル の小型化が課題となっている。また、情報の「0」「1」を切替える方法にも、磁界やスピントル ク注入による磁化反転、磁壁移動があるが、情報書き込み速度の高速化も課題の一つである。 近年、これらの課題を解決するために「ジャロシンスキー・守谷相互作用 (DMI)」[1, 2] が注目 されている。DMI は磁化構造を捻らせる作用を持ち、スカーミオン構造 [3] や高速磁壁移動 [4] といった特異な磁性現象が生じる。スカーミオン構造は、磁気渦構造体のひとつで、直径が 10 nm程度と極めて小さい構造であるため、記憶セルの小型化が期待されている。また磁壁移動 により情報書き込みを行う「磁壁移動型 MRAM」において、高速磁壁移動は情報書き込み速 度の高速化が期待できる。そのため DMI および DMI による磁性現象は近年盛んに研究が行わ れている。DMI の強度 (D:DMI 定数) を測定する手法についても研究が進められ、DMI 定数の 測定には直接磁性体の磁化の捻れを測定する方法 [5, 6] や、DMI による磁性現象を利用した測 定法 [4, 7] などがある。
1.2
目的
本研究では、DMI 定数の測定法として新たに磁化反転実験により測定する手法を提案し、マ イクロマグネティックシミュレーションにて検討した。具体的には反転磁界および反転確率を 用いる 2 つの手法による測定法を提案し、シミュレーションにより DMI による磁性現象への影 響を調査、解析し、提案手法の有効性を確認した。1.3
論文構成
本論文は、以下の章で構成されている。 第 2 章 本研究での基本事項について述べる。 第 3 章 マイクロマグネティックシミュレーションで用いる数値解法について述べる。 第 4 章 本研究で用いるシミュレーションモデルについて述べる。 第 5 章 本研究で検討する測定法について述べる。 第 6 章 反転磁界による測定法でのシミュレーション結果について述べる。 第 7 章 反転確率の外部磁界依存性による測定法でのシミュレーション結果について述べる。 第 8 章 反転確率の外部磁界印加角度によるシミュレーション結果について述べる。 第 9 章 本論文のまとめについて述べる。第
2
章 基本事項
第 2 章では、本研究での基本となる用語について説明をする。まず本研究で行うマイクロマ グネティックシミュレーション、磁性体を取り巻く物理的原理及び熱エネルギーによる磁性体 への影響について説明する。そして本研究の主題であるジャロシンスキー・守谷相互作用及び 本作用の従来の測定法について説明する。2.1
マイクロマグネティックシミュレーション
磁石を分割していくと図 2.1 のように原子レベルまで分割できるとするモデルが考えられる。 このモデルの場合、各原子は N 極と S 極を持つ磁気モーメント (後述) を持つと考えられる。これ らの磁気モーメントによって作られる磁性体内部の磁化構造を扱う分野をマイクロマグネティッ クスといい、これらの現象をシミュレーションにより解析することをマイクロマグネティック シミュレーションという。マイクロマグネティックシミュレーションの特徴として、実験と比 較して短時間かつ低コストで実験を行うことができるという点や、実験的には困難であるナノ 秒単位での磁化構造の変化を詳細に観測できるというメリットがある。 図 2.1: 永久磁石の分解図2.2
磁気モーメント
磁気モーメントとは、磁極の対を表す物理量で、単位は emu である。原子が持つ磁気モーメ ントのことを原子磁気モーメントと呼び、単位体積当たりの量を磁化と呼ぶ。 図 2.2 のように 磁極間の距離を l、磁極の強さを q とすると、磁気モーメントの大きさは ql となる。図 2.2: 磁気モーメントの概念図
2.3
LLG
方程式
LLG(Landau-Lifshitz-Gilbert)方程式とは、磁界による原子磁気モーメントの運動を表す方 程式である。LLG 方程式は、式 (2.1) で表される [8]。 ˙~ M = −|γ|( ~M × ~H) + α M( ~M ×M )˙~ (2.1) γ :磁気回転比 α :損失定数 ~ M :原子磁気モーメント ~ H :磁界 式 (2.1) において第 1 項をジャイロ項と呼び、磁界周りの磁気モーメントの歳差運動を表し、 第 2 項を損失項と呼び、磁気モーメントが磁界方向に緩和する運動を表す。図 (2.3) に歳差運動 と緩和運動による原子磁気モーメントの運動を示す。 図 2.3: 原子磁気モーメントの運動2.4
ジャロシンスキー・守谷相互作用
(DMI)
ジャロシンスキー・守谷相互作用 (DMI)[1, 2] とは、隣接する原子磁気モーメントがねじれる 作用である。DMI により生じるエネルギーは式 (2.2) で表される [9]。 EDM I = − ~Dij · ~Si × ~Sj (2.2) ~ Dij : DMIベクトル、~Si, ~Sj :隣接する磁気モーメント DMIは、 「バルクの非対称磁性材料」および「強力なスピン軌道結合を有する強磁性体と 非磁性金属との間の界面」で観測されている。これらの DMI はそれぞれ「バルク DMI[10]」、 「界面 DMI[11]」と呼ばれ、性質が異なる。DMI を用いた磁性現象には、スカーミオン [3] や磁 壁移動の高速化 [4] がある。スカーミオン構造は、最初はバルクの非対称磁性材料で観測された が近年超薄膜で観測されている。これは超薄膜では強力なスピン軌道結合を有する隣接層への 近接によって誘発される界面 DMI によって説明されている [9]。また界面 DMI は強磁性体や非 磁性金属の膜厚を変えることにより DMI の強度が実験的に可変 [12] であるとされている。本研 究では界面 DMI に着目する。2.4.1
界面ジャロシンスキー・守谷相互作用
(
界面
DMI)
前述のように、界面ジャロシンスキー・守谷相互作用 (界面 DMI) とは、「強力なスピン軌道 結合を有する強磁性体と非磁性金属との間の界面」にて観測される DMI である。「強力なスピ ン軌道結合を有する強磁性体」の例に Co[12],Fe[9],CoFeB[13] がある。また「非磁性金属」の例 として、Pt[12],W[13],Ir[9] などがある。2.4.2
DMI
の測定
DMI定数の定量的評価は種々の方法がある。測定法は、観測される磁化構造により直接測定 する方法 (X 線回折 [5]、ブルリアン散乱 [6]) や、DMI 定数により変化する磁性現象 (スカーミ オンの半径 [7]、磁壁移動速度 [4] など) で推定する方法がある。2.5
実効磁界
実効磁界とは、原子磁気モーメントが受ける各磁気エネルギーを変分することにより磁界に 変換したものである。磁気エネルギーを変分し実効磁界を求める式は、(2.3) で表される。本研 究で用いる実効磁界は、異方性磁界、交換磁界、静磁界、DMI 磁界、外部磁界で構成されるも のとした。異方性磁界、交換磁界、静磁界、DMI 磁界については以下の項目で説明する。 ~ H = − δǫ δ ~M (2.3) ~ H :実効磁界, ǫ : 磁気エネルギー, ~M :原子磁気モーメント2.5.1
異方性磁界
異方性磁界とは、原子磁気モーメントがある特定方向に向こうとする性質による磁気エネル ギーを磁界に変換したものである。原子磁気モーメントが一つの軸方向に向こうとする性質を 一軸磁気異方性といい、向こうとする一つの軸を磁化容易軸という。磁化容易軸が z 軸である 場合の磁気エネルギーは式 (2.4) で表される。ここで Kuは磁気異方性定数 (単位は erg/cm3)で ある。 ǫK = Ku(1 − m2z) (2.4)2.5.2
交換磁界
交換磁界とは、隣接する原子磁気モーメント間の角度差によって生じる磁気エネルギーを磁 界に変換したものである。ここでの磁気エネルギーを交換エネルギーという。交換エネルギー は、式 (2.5) で表される。式中の A は交換スティッフネス定数 (単位は erg/cm) である。 ǫA= A(∇~m)2 (2.5)2.5.3
静磁界
静磁界とは、磁性体が作り出す磁界である。静磁界は、磁性体内の原子磁気モーメントの配 列によって決まる。一般的に、磁気モーメントの配列は解析的に表すことができないため、一 般的に静磁界も解析的に表すことができない。2.5.4
DMI
磁界
DMI磁界とは、DMI によって生じる磁気エネルギーを磁界に変換したものである。2.6
熱揺らぎ
熱揺らぎとは、熱エネルギーにより原子磁気モーメントが振動する現象である。熱揺らぎは ~ HT として式 (2.6) で表される [14]。 hHiT(t)H T j (t + τ )i = 2kBT α |γ|vMs δ(τ )δij (2.6) kB :ボルツマン定数, T : 温度 (K), α : 損失定数 γ :磁気回転比 (rad/(s · Oe)), v : 計算セルの体積 (cm3) Ms :飽和磁化 (emu/cm3), 熱揺らぎの効果を考慮した LLG 方程式は、式 (2.7) で表される [14]。 ˙~ M = −|γ| ~M × ~H + ~HT+ α M ~M ×M˙~ (2.7)第
3
章 マイクロマグネティックシミュレー
ション
第 3 章では、マイクロマグネティックシミュレーションで用いる数値解法について説明する。 まず LLG 方程式及び実効磁界の数値解法について説明する。次に、実効磁界計算の高速化手法 や、熱揺らぎ磁界の乱数生成方法について説明する。3.1
LLG
方程式の数値解法
LLG方程式は式 (3.1) で表される (再掲)。 ˙~ M = −|γ|( ~M × ~H) + α M( ~M ×M )˙~ (3.1) γ :磁気回転比 α :損失定数 ~ M :原子磁気モーメント ~ H :磁界 ここで式 (3.1) は両辺に ˙~Mが含まれているため、このままの形ではシミュレーションで使用 することができない。そのためにシミュレーションで使用できる形に変形する必要がある。ま ず ~M = M ~mとおき、式 (3.1) を整理すると、式 (3.3) となる。 M ˙~m = −|γ|M ~m × ~H + α M(M 2~ m × ˙~m) (3.2) ⇐⇒ ˙~m = −|γ|~m × ~H + α( ~m × ˙~m) (3.3) 次に、両辺に ~mの外積を取り、整理すると式 (3.4) となる。 ~ m × ˙~m = −|γ|~m × (~m × ~H) + α ~m × (~m × ˙~m) = −|γ|[( ~H · ~m) ~m − (~m · ~m) ~H] + α[( ˙~m · m)~m − (~m · ~m) ˙~m] = −|γ|[( ~H · ~m) ~m − ~H] + α(−m) (3.4) さらに、式 (3.3) の右辺中の ~m × ˙~mに式 (3.4) を代入して、整理すると式 (3.7)˙~m = −|γ|~m × ~H + α[−|γ|(( ~H · ~m)~m − ~H) + α(− ˙~m)] = −|γ|~m × ~H − α|γ|[( ~H · ~m) ~m − ~H] − α2 ˙~m (3.5) ⇐⇒ (1 + α2) ˙~m = −|γ|~m × ~H − α|γ|[( ~H · ~m) ~m − ~H] (3.6) ⇐⇒ ˙~m = −|γ| 1 + α2[ ~m × ~H − α(( ~H · ~m) ~m − ~H)] (3.7) また d = −|γ| 1+α2, e = ~H · ~mとおくと式 (3.7) は式 (3.8) で表される。 ˙~m = d(~m × ~H + α(e~m − ~H)) (3.8) 式 (3.8) を x, y, z 成分毎で表すと、式 (3.9) ∼ (3.11) のようになる。 ˙ mx= d((myHz− mzHy) + α(emx− Hx)) (3.9) ˙ my = d((mzHx− mxHz) + α(emy− Hy)) (3.10) ˙ mz = d((mxHy − myHx) + α(emz− Hz)) (3.11) シミュレーションでは、式 (3.9),(3.10),(3.11) を利用する。
3.2
オイラー法
オイラー法とは、LLG 方程式などの常微分方程式の解を求めるための手法である。まず、オ イラー法による微分方程式の数値解法について式 (3.12) を用いて説明する。 dy dt = f (t, y(t)) (3.12) この式を数値的に解くために、時間を ∆t で離散化する。y(t) の 1 ステップ後の値 y(t + ∆t) は Taylor 展開より式 (3.13) で表される。 y(t + ∆t) = y(t) + y′(t)∆t +1 2y ′′(t)∆t2+ 1 6y ′′′(t)∆t3+ O(∆t4) (3.13) ここで ∆t は極めて小さいものとし、∆t2以降の項を打ち切ると、式 (3.14) が得られる。 y(t + ∆t) = y(t) + y′(t)∆t (3.14)式 (3.14) に従い、y の時間変化を追いかける手法をオイラー法と呼ぶ。オイラー法を式 (3.9) ∼ (3.11) に対して適用すると、式 (3.15),(3.16),(3.17) となる。 mx(t + ∆t) = mx(t) + ∆t · d((myHz− mzHy) + α(emx− Hx)) (3.15) my(t + ∆t) = my(t) + ∆t · d((mzHx− mxHz) + α(emy− Hy)) (3.16) mz(t + ∆t) = mz(t) + ∆t · d((mxHy− myHx) + α(emz− Hz)) (3.17) オイラー法では 1 回の計算による打ち切り誤差が O(∆t2)となる。そのため、一定時間のシ ミュレーションを行った場合の大域誤差は O(∆t) となる。このため、オイラー法は一次の手法 と呼ばれる。
3.3
4
次のルンゲ・クッタ法
3.2で述べたオイラー法は一次の精度であった。本研究では、より高次の手法である 4 次のル ンゲ・クッタ法を用いた。4 次のルンゲ・クッタ法は式 (3.18) で表される。ここで式 (3.18) 中 の k1, k2, k3, k4は式 (3.19) で表される。 ~ m(t + ∆t) = ~m(t) + 1 6(k1 + 2k2+ 2k3+ k4) (3.18) k1 = f (t, ~m(t)) k2 = f (t + ∆t2 , ~m(t) + k21) k3 = f (t + ∆t2 , ~m(t) + k22) k4 = f (t + ∆t, ~m(t) + k3) (3.19) ここで式 (3.19) 中の f(t, ~m(t))は式 (3.20) で表される。 f (t, ~m(t)) = ∆t · (~m × ~H + α(d ~m(t) − ~H)) (3.20) 式 (3.18) の計算で生じる打ち切り誤差は O(∆t5)となるので、大域誤差は O(∆t4)となる [15]。3.4
計算領域の離散化
本研究では計算領域を二次元モデルで離散化した。二次元モデルでの計算領域の離散化の概 要図を図 3.4 で表す。図中の直方体を計算セルとし、計算セル一つに対して原子磁気モーメン トを一つ配置する。図 3.1: 計算領域の離散化の概要図
3.5
実効磁界の計算
本研究における実効磁界は異方性磁界、交換磁界、DMI 磁界、静磁界、外部磁界からなるも のとした (式 (3.21))。 ~ H = ~HK+ ~HA+ ~HDMI+ ~HD + ~HEXT (3.21) HK :異方性磁界 HA:交換磁界 HDMI: DMI磁界 HD :静磁界 HEXT :外部磁界3.5.1
異方性磁界
異方性エネルギーは、式 (3.22) で表される (再掲)。 ǫK = Ku(1 − m2z) (3.22) 式 (3.22) を変分することにより異方性磁界が求まる。 ~ HK = −δǫ K δ ~M = − 1 Ms δǫK δ ~m (3.23) 式 (3.23) に式 (3.22) を代入し、各成分毎に分解し、整理すると式 (3.24) ∼ (3.26) となる。 HK x = 1 Ms − ∂ǫK ∂mx + ∂ ∂x ∂ǫK ∂ ∂mx ∂x ! + ∂ ∂y ∂ǫK ∂∂mx ∂y = 0 (3.24) HyK = 1 Ms − ∂ǫK ∂my + ∂ ∂x ∂ǫK ∂∂my ∂x + ∂ ∂y ∂ǫK ∂∂my ∂y = 0 (3.25)HK z = 1 Ms − ∂ǫK ∂mz + ∂ ∂x ∂ǫK ∂ ∂mz ∂x ! + ∂ ∂y ∂ǫK ∂∂mz ∂y = 2Ku Ms mz (3.26) 式 (3.24) ∼ (3.26) より、異方性磁界は式 (3.27) で表される。 ~ HK = 0 0 2Ku Ms mz (3.27)
3.5.2
交換磁界
交換エネルギーは式 (3.28) で表される (再掲)。 ǫA= A(∇~m)2 (3.28) ここで二次元モデルの場合、交換エネルギーは式 (3.29) で表される。 ǫA= A ( ∂mx ∂x 2 + ∂my ∂x 2 + ∂mz ∂x 2 + ∂mx ∂y 2 + ∂my ∂y 2 + ∂mz ∂y 2) (3.29) 式 (3.29) を変分することにより交換磁界が求まる。 ~ HA = −δǫ A δ ~M = − 1 Ms δǫA δ ~m (3.30) 式 (3.30) に式 (3.29) を代入し、各成分毎に分解すると、式 (3.31) ∼ (3.33) となる。 HxA = 1 Ms − ∂ǫA ∂mx + ∂ ∂x ∂ǫA ∂ ∂mx ∂x ! + ∂ ∂y ∂ǫA ∂∂mx ∂y = 2A Ms ∂2m x ∂x2 + ∂2m x ∂y2 (3.31) HA y = 1 Ms − ∂ǫA ∂my + ∂ ∂x ∂ǫA ∂∂my ∂x + ∂ ∂y ∂ǫA ∂∂my ∂y = 2A Ms ∂2m y ∂x2 + ∂2m y ∂y2 (3.32)HA z = 1 Ms − ∂ǫA ∂mz + ∂ ∂x ∂ǫA ∂ ∂mz ∂x ! + ∂ ∂y ∂ǫA ∂∂mz ∂y = 2A Ms ∂2m z ∂x2 + ∂2m z ∂y2 (3.33) ここで式 (3.31) ∼ (3.33) の中にある二階偏微分は、中心差分法を用いて計算した。 ∂2m x ∂x2 ≃
mxi+1,j − 2mxi,j + mxi−1,j
∆x2 (3.34)
∂2m y
∂x2 ≃
myi+1,j − 2myi,j + myi−1,j
∆x2 (3.35)
∂2m z
∂x2 ≃
mzi+1,j − 2mzi,j + mzi−1,j
∆x2 (3.36)
∂2m x
∂y2 ≃
mxi,j+1 − 2mxi,j + mxi,j−1
∆y2 (3.37)
∂2m y
∂y2 ≃
myi,j+1 − 2myi,j + myi,j−1
∆y2 (3.38)
∂2m z
∂y2 ≃
mzi,j+1 − 2mzi,j + mzi,j−1
∆y2 (3.39) 式 (3.34) ∼ (3.39) より、交換磁界は式 (3.40) で表される。 ~ HA= 2A Ms 1 ∆x2
mxi+1,j− 2mxi,j + mxi−1,j
myi+1,j − 2myi,j + myi−1,j
mzi+1,j − 2mzi,j + mzi−1,j
+ 1 ∆y2
mxi,j+1 − 2mxi,j + mxi,j−1
myi,j+1 − 2myi,j + myi,j−1
mzi,j+1− 2mzi,j + mzi,j−1
(3.40)
3.5.3
DMI
磁界
一般的に隣接する磁気モーメント Si, Sj間で働く DMI エネルギーは式 (3.41) で表される (再 掲)[9]。 EDMI = − ~Dij · ~Si× ~Sj (3.41) 薄膜の場合、DMI エネルギーは式 (3.42) で表される [9] ǫDMI= −D mx ∂mz ∂x − mz ∂mx ∂x + my ∂mz ∂y − mz ∂my ∂y (3.42)式 (3.42) 中の D は DMI 定数であり、DMI ベクトル ~Dijの大きさに対して 1/at に比例した面積 あたりのエネルギー密度である (a は格子定数, t は薄膜の厚さ)[9]。式 (3.42) を変分することで DMI磁界が求まる。 ~ HDMI = −δǫ DMI δ ~M = − 1 Ms δǫDMI δ ~m (3.43) 式 (3.43) に式 (3.42) を代入し、各成分毎に分解すると、式 (3.44) ∼ (3.46) が得られる。 HxDMI = 1 Ms − ∂ǫDMI ∂mx + ∂ ∂x ∂ǫDMI ∂ ∂mx ∂x ! + ∂ ∂y ∂ǫDMI ∂∂mx ∂y + ∂ ∂z ∂ǫDMI ∂ ∂mx ∂z ! = 2D Ms ∂mz ∂x (3.44) HyDMI = 1 Ms − ∂ǫDMI ∂my + ∂ ∂x ∂ǫDMI ∂∂my ∂x + ∂ ∂y ∂ǫDMI ∂∂my ∂y + ∂ ∂z ∂ǫDMI ∂∂my ∂z = 2D Ms ∂mz ∂y (3.45) HzDMI = 1 Ms − ∂ǫDMI ∂mz + ∂ ∂x ∂ǫDMI ∂ ∂mz ∂x ! + ∂ ∂y ∂ǫDMI ∂∂mz ∂y + ∂ ∂z ∂ǫDMI ∂ ∂mz ∂z ! = −2D Ms ∂mx ∂x + ∂my ∂y (3.46) ここで、式 (3.44) ∼ (3.46) の中にある一階偏微分は、中心差分を用いて計算した。 ∂mx ∂x ≃ mxi+1,j − mxi−1,j 2∆x (3.47) ∂my ∂y ≃ mxi,j+1 − mxi,j−1 2∆y (3.48) ∂mz ∂x ≃ mxi+1,j − mxi−1,j 2∆x (3.49) ∂mz ∂y ≃ mxi,j+1 − mxi,j−1 2∆y (3.50) 式 (3.47) ∼ (3.50) より、DMI 磁界は式 (3.51) で表される。 HDM I = 2D Ms −mz i+1−mzi −1 2∆x −mzj+12∆y−mzj−1 mxi+1−mxi−1 2∆x + myj+1−myj −1 2∆y (3.51)
3.5.4
境界条件
DMIは、隣接する磁気モーメントを捻らせる作用を持つ。そのため、本研究では、DMI に より境界部が捻れることを考慮した境界条件を用いた。DMI を考慮した境界条件の式は、式 (3.52)で表される [16]。 ∂ ~m ∂~n = D 2A(ˆz × ~n) × ~m (3.52) 式 (3.52) 中の ~n は適用する境界部の方向ベクトルを表す。例えば x 方向の境界で適用する場 合の ~n は ˆx である。式中の一階偏微分は、前進差分近似を用いて表すとする (式 (3.53))。 ∂ ~m ∂~n = ~ mi+1,j− ~mi,j ∆n (3.53) ここで各方向の境界で用いる境界条件を求める。まず x 方向の境界での境界条件を求めるた めに、式 (3.52) 中の ~n に ˆx を代入すると、式 (3.54) が得られる。 ∂ ~m ∂ ˆx = D 2A(ˆz × ˆx) × ~m = D 2A mz 0 −mx (3.54) 式 (3.53) と式 (3.54) より、境界部の外にある仮想的な計算点 (i = −1 or Nx) での磁気モーメン トの向きは、式 (3.55) ∼ (3.60) で表される。mxi=−1,j = mxi=0,j− ∆x · mzi=0,j ·
D
2A (3.55)
myi=−1,j = myi=0,j (3.56)
mzi=−1,j = mzi=0,j + ∆x · mxi=0,j ·
D
2A (3.57)
mxi=Nx,j = mxi=Nx−1,j + ∆x · mzi=Nx−1,j ·
D
2A (3.58)
myi=Nx,j = myi=Nx−1,j (3.59)
mzi=Nx,j = mzi=Nx−1,j − ∆x · mxi=Nx−1,j ·
D 2A (3.60) 次に、y 方向の境界での境界条件を求めるために、式 (3.52) 中の ~n に ˆyを代入すると、式 (3.61) が得られる。 ∂ ~m ∂ ˆy = D 2A(ˆz × ˆy) × ~m = D 2A 0 mz −my (3.61)
式 (3.53) と式 (3.61) より、境界部の外にある仮想的な計算点 (j = −1 or Ny) での磁気モー メントの向きは、式 (3.62) ∼ (3.67) で表される。
mxi,j=−1 = mxi,j=0 (3.62)
myi,j=−1 = myi,j=0 − ∆y · mzi,j=0 ·
D
2A (3.63)
mzi=,j=−1 = mzi,j=0 + ∆y · myi,j=0 ·
D
2A (3.64)
mxi,j=Ny = mxi,j=Ny −1 (3.65)
myi,j=Ny = myi,j=Ny −1 + ∆y · mzi,j=Ny −1 ·
D
2A (3.66)
mzi,j=Ny = mzi,j=Ny −1 − ∆y · myi,j=Ny −1 ·
D 2A (3.67)
3.5.5
静磁界
静磁界とは、磁性体が作り出す磁界である。 磁性体内部の磁化構造は磁性体内部での相互作 用や外部から受ける磁界などにより変化するため、静磁界を解析的に求めることは不可能であ り、数値的に求めることにする。計算領域は 3.4 で述べたように、直方体セルを用いて離散化 を行う。 一つの計算セル内では磁気モーメントは全て同じ方向を向くと仮定すると、これによって計 算領域内に現れる磁極は全て計算セルの表面に現れることになる。一つの計算点での静磁界は、 計算領域内のすべての磁極が今考えている計算点に作り出す静磁界の和として求めることがで きる。 計算セルの各面は xy,yz,xz 面にそれぞれ平行であるとする。計算セルの各面は、それぞれ面 密度 ±Mx, ±Myおよび ± Mzで磁極は分布しているものとする。これらの 6 つの面を、yz 面に 平行な面、xz 面に平行な面および xz 面に平行な面を分け、それぞれの面に現れる磁極が観測 点に作る静磁界を求める。 まず yz 面に平行な右側の面を考える。観測点 (x1, y, z)離れた面上の微小領域が観測点に作 り出す磁界は式 (3.68) ∼ (3.70) で表される。なお、r =px2 1+ y2+ z2とする。 ∆HD x = − Mx r2 x1 r ∆y∆z (3.68) ∆HD y = − My r2 y r∆y∆z (3.69) ∆HD z = − Mz r2 z r∆y∆z (3.70)図 3.2: 計算領域の離散化 対象とする面上の磁極が観測点に作り出す磁界は、これらの磁界を扱う面 (x1, y, z)にわたって 式 (3.68),(3.69),(3.70) を積分することによって求められる。 ∆HxD = Z z1 z0 Z y1 y0 ∆HxD (3.71) ∆HD y = Z z1 z0 Z y1 y0 ∆HD y (3.72) ∆HzD = Z z1 z0 Z y1 y0 ∆HzD (3.73) 同様に yz 面に平行な左側の面上の磁極が観測点に作り出す磁界を求め、これらの式をま とめる。同様の操作を xz 面に平行な面で行い、それぞれをまとめる。計算セル上での磁荷 が観測点 (i, j) に作り出す磁界は式 (3.74) ∼ (3.76) で表される。式 (3.74),(3.75),(3.76) 中の qxx, qxy, qxz, qyy, qyz, qzzは静磁界係数と呼ぶ。
HD x (i, j) = nx X i0=1 ny X j0=1 [qxx(i0, j0) · mx(i0− i, j0− j) +qxz(i0, j0) · my(i0− i, j0− j) +qxz(i0, j0) · mz(i0− i, j0− j)] (3.74) HyD(i, j) = nx X i0=1 ny X j0=1 [qxy(i0, j0) · mx(i0− i, j0− j) +qyy(i0, j0) · my(i0− i, j0− j) +qyz(i0, j0) · mz(i0− i, j0− j)] (3.75)
HD z (i, j) = nx X i0=1 ny X j0=1 [qxz(i0, j0) · mx(i0− i, j0− j) +qyz(i0, j0) · my(i0− i, j0− j) +qzz(i0, j0) · mz(i0− i, j0− j)] (3.76) 式 (3.74),(3.75),(3.76) は convolution 演算と呼ばれ、nx, ny はそれぞれ x, y 方向の計算点の数 である。ここまでの静磁界計算では、静磁界が計算セル中央部に現れるものと仮定して計算し た。しかし実際には計算セル内でも静磁界が変化している。そのため、より正確な静磁界を求 めるためには、計算セル内に現われる静磁界の平均値を用いる必要がある。計算セル内で平均 した静磁界を求めるための静磁界係数を以下に示す [17]。 F 1(X; Y, Z) = XY Z arctan Y Z XD +1 2Y (Z 2 − X2) ln |D − Y | +1 2Z(Y 2 − X2) ln |D − Z| + 1 6(Y 2+ Z2 − 2X2)D (3.77) F 2(Z; X, Y ) = −XY Z ln |D + Z| + 16Y (Y2− 3Z2) ln |D + X| +1 6X(X 2 − 3Z2) ln |D + Y | + 12X2Z arctan Y Z XD +1 2Y 2Z arctan XZ Y D +1 6Z 3arctan XY ZD +1 3XY D (3.78) < qxx >v= 1 v 3 X i,j,k=1
(−1)i+j+k−1sn(i)sn(j)sn(k) × F 1 (x + ax(i); y + ay(j), z + az(k))(3.79)
< qxy >v= 1 v 3 X i,j,k=1 (−1)i+j+k−1
sn(i)sn(j)sn(k) × F 2 (z + az(k); x + ax(i), y + ay(j))(3.80)
< qyy >v= 1 v 3 X i,j,k=1 (−1)i+j+k−1
sn(i)sn(j)sn(k) × F 1 (y + ay(j); z + az(k), x + ax(i))(3.81)
< qzz >v= 1 v 3 X i,j,k=1 (−1)i+j+k−1
sn(i)sn(j)sn(k) × F 1 (z + az(k); x + ax(i), y + ay(j))(3.82)
< qxz >v= 1 v 3 X i,j,k=1 (−1)i+j+k−1
< qyz >v= 1 v 3 X i,j,k=1
(−1)i+j+k−1sn(i)sn(j)sn(k) × F 2 (x + ax(i); y + ay(j), z + az(k))(3.84)
ここで D, ax(1), ax(2), ax(3), ay(1), ay(2), ay(3), az(1), az(2), az(3), sn(1), sn(2), sn(3) は式 (3.85),(3.86) で表される。
D =√X2+ Y2+ Z2 (3.85)
ax(1) = −ddx, ax(2) = 0, ax(3) = ddx ay(1) = −ddy, ay(2) = 0, ay(3) = ddy az(1) = −ddz, az(2) = 0, az(3) = ddz
sn(1) = 1, sn(2) = 2, sn(3) = 1 (3.86) 式 (3.86) 中の、ddx, ddy, ddz は計算セルの各座標軸における長さである。 式 (3.74),(3.75),(3.76) の計算量は、計算点 n に対して O(n2)となる。そのため静磁界計算の 計算量は O(n2)となり、LLG 方程式の数値計算において計算時間の大部分が静磁界に計算に費 される。3.6 では、静磁界計算の高速化について述べる。
3.6
静磁界の高速計算
(
離散高速フーリエ変換
)
前述のように、LLG 方程式の数値計算では、計算量が O(n2)となる静磁界計算に計算時間の大 部分が費され、このままでは計算点が一万点以上の大規模な計算を行うことは実用上不可能であ る。そこで離散高速フーリエ変化を用いることによって、静磁界の計算式に現れる Convolution 演算を高速に行う手法について述べる。3.6.1
離散高速フーリエ変換を用いた
Convolution
演算の高速計算
式 (3.87) について考える。 HD(i) = n X j=1 Cq(j − i) · M(j), (i = 1, · · · , n) (3.87) ここで、HD, M はそれぞれ長さ n のベクトル、Cq は長さ 2n − 1 のベクトル、n は二の累乗の 数とする。この演算は Convoluton 演算と呼ばれており、HD (1), · · · , HD(n)を全て求めるため に必要な計算時間は n の二乗に比例する。次に M および Cq が周期 n の周期関数である場合を 考える。周期性より Cq は長さ n のベクトルで表すことができる。このとき、上記の演算は離 散フーリエ変換を用いて、以下の手順で計算することができる。[18, 19] 1. Mおよび Cq のフーリエ成分 M∗, Cq∗を求める。 2. 求めたフーリエ成分を掛け合わせ、結果を HD∗とした。3. HD∗を逆フーリエ変換する。このとき得られたものの実数部が求めるべき HDである。 上記の手順に従い計算を行ったときの計算時間を説明する。手順 1 及び 3 の長さ n のベクトル の離散フーリエ変換および逆変換であるため、離散高速フーリエ変換の手法を用いて計算を行 えば、n log(n) に比例する計算時間で計算を行うことができる。手順 2 は n に比例する計算時 間で計算を行うことができる。従って、全ての計算は O(n log(n)) で行うことができると考えら れる。図 3.6.1 中の直接法は 3.5.5 で示した数値解法のことである。 図 3.3: 直接法と FFT の計算量の関係
3.6.2
非周期構造での
Convolution
演算
上記の計算は計算対象が周期構造を持つ場合に適応することができる手法であるが、計算対 象が周期構造を持たない一般の場合に対しても、以下に述べるゼロパディング手法を用いるこ とにより適用することができる [19]。 まず長さ n のベクトル M および Cq を以下のように長さ 2n のベクトルに拡張し、それぞれ M′, Cq′とする。 M′ : M (1), M (2), · · · , M(n), 0, 0, · · · , 0 Cq′ : 0, Cq(−n + 1), Cq(−n + 2), · · · , Cq(0), Cq(1), · · · , Cq(n − 1) ベクトル M’ は、M の成分に対して、要素が 0 の n 個の成分を付け加えたものとする。ベク トル Cq’ は、元々の Convolution 演算に必要な 2n − 1 個の要素と、ひとつの 0 の成分から成る ものとする。このようにして得られた M′と Cq′に対して、前節の手順に従い演算を行う。得ら れたベクトルを HD’とすると、HD ’(n + 1), · · · , HD’(2n)の成分が求める演算の解である。上 記の例では、ベクトル M’ を作成するときに、後半の n 個の成分を 0 としているが、逆に前半 の n 個の成分を 0 にすることも可能である。この場合求める HDは、得られた HD’の前半の成 分となる。 要素数 n が二の累乗でない場合は、n より大きい二の累乗の中で最小のものを求め、この数 を新に n として計算を行う。ここでベクトル M の拡張された領域の要素は 0 とする。図 3.4: 2 次元モデルでのフーリエ変換の概要図
3.7
熱揺らぎ磁界における乱数生成
熱揺らぎの効果は ~HT として実効磁界に換算することで LLG 方程式に組み込むことが可能で ある。 ~HT は式 (3.88) で表される (再掲)。 hHT i (t)HjT(t + τ )i = 2kBT α |γ|vMs σ(τ )σij (3.88) 式内の hi は時間平均を表す。HT は平均値が 0 で標準偏差が σ となる正規分布を持つ乱数で 表す。 σ = s 2kBT γvMs∆t (3.89) kB :ボルツマン定数 T :温度 (K) γ :磁気回転比 (rad/(s · Oe)) v :計算セルの体積 (cm3) Ms:飽和磁化 (emu/cm3) ∆t :時間刻み (s) この乱数は、MT 法 [20] により一様乱数を生成し、それを Box-muller 法により正規分布となる 乱数に変換し、生成した。~ HT を LLG 方程式に組み込むと式 (3.90) となる [14]。 ˙~ M = −|γ| ~M × ~H + ~HT+ α M ~M ×M˙~ (3.90)
第
4
章 シミュレーションモデル
第 4 章では、本研究で用いたシミュレーションモデルについて説明する。まずマクロスピン モデルと、マクロスピンモデルで用いた材料定数について説明し、マイクロマグネティックモ デルにおける離散化方法と、用いた材料定数について説明する。4.1
マクロスピンモデル
マクロスピンモデルとは、磁性体内部の磁化構造を一様と仮定したモデルである。そのため 交換磁界及び DMI 磁界といった磁気モーメント間の相互作用は生じない。また磁性体内部の磁 化構造を一様と仮定していることから、シミュレーション結果と解析式の結果を直接比較可能 なモデルである。本研究では熱揺らぎを考慮した磁化反転シミュレータの妥当性の検討のため にマクロスピンモデルを用いた (後述)。マクロスピンモデルにおける実効磁界は異方性磁界と 外部磁界で構成され、式 (4.1) で表される。 図 4.1: マクロスピンモデル ~ H = ~HK+ ~HEXT (4.1) マクロスピンモデルで用いた材料定数はを以下に示す。想定材料は CoFeB である [21]。 飽和磁化 : Ms= 1500 emu/cm3 磁気異方性定数 : Ku = 14.0 Merg/cm3 磁気回転比 : γ = 1.76 × 107 rad/(s · Oe) 損失定数 : α = 1.0 また、熱揺らぎを考慮したマクロスピンモデルの場合、熱エネルギーの影響により磁気モー メントが z 軸に対して θ0(rad)傾いた状態が平衡状態となることが知られている。[23]。θ0 = r 1 ∆ (rad) (4.2) ∆ = KuV kBT (4.3) 式 (4.3) 中の ∆ は磁気モーメントの熱に対する耐久性を表す熱安定性指数である。熱揺らぎを 考慮した磁化反転シミュレータの妥当性を検討するために、平衡状態の磁気モーメントの平均 磁化角度を式 (4.2) と比較する (後述)。
4.2
マイクロマグネティックモデル
マイクロマグネティックモデルとは、磁性体内部の磁化構造を複数の原子磁気モーメントを 用いて表現するモデルで、実際の磁性体に近いモデルである。本研究で用いたマイクロマグネ ティックモデルでは円盤ディスクを直方体セルで離散化したものを用いた。セル 1 つのサイズ は、2 nm × 2 nm × 1 nm とした。 図 4.2: 円盤の離散化モデル マイクロマグネティックモデルで用いた材料定数を以下に示す。想定材料はマクロスピンモ デル同様、CoFeB[21] である。 飽和磁化 : Ms = 1500 emu/cm3 磁気異方性定数 : Ku = 14.0 Merg/cm3 磁気回転比 : γ = 1.76 × 107 rad/(s · Oe) 交換スティフネス定数 : A = 3.1 µerg/cm 損失定数 : 後述 DMI定数 : 後述第
5
章 提案する
DMI
定数の測定法
第 5 章では、外部磁界を用いて DMI 定数を測定する手法を 2 つ提案する。最初の手法は反 転磁界を用いる手法であり、次の手法は反転確率を用いる手法である。以下に各手法の説明を する。5.1
反転磁界による測定
まず反転磁界による測定について述べる。DMI により磁性体内部の磁化構造が変化する。ま た円盤直径により DMI 定数による磁化構造の変化量が変わることが考えられる。これらの磁 化構造の変化により、反転磁界が変化することが考えられる。また反転磁界の変化量は、面内 反転磁界と面直反転磁界で異なることが考えられる。これらのことより、直径が異なる磁性円 盤の面内及び面直反転磁界より、DMI 定数を測定できる可能性が考えられる。この手法による 調査手順を以下に示す。まず初期磁化を面直上向きの飽和状態と設定する。これに、面直下向 きと面内方向にそれぞれ外部磁界 (HEXT z , HxEXT)を印加することにより面直反転磁界 (Hzsw)と 面内反転磁界 (Hsw x )を円盤直径と DMI を変えて計算する。得られた結果より DMI 定数を推定 する。 図 5.1: 反転磁界による測定法の概略図5.2
反転確率による測定
次に反転確率による測定手法について述べる。熱揺らぎを考慮した場合、熱揺らぎのランダ ム性のために原子磁気モーメントの磁化反転は確率的に起こる [22]。このため反転磁界ではな く、反転確率を求めることになる。DMI により磁性体内部の磁化構造が変化するため、磁界に よる反転確率が変化する可能性が考えられる。磁界パルス幅を変えることにより、磁界による反転確率が変化することも考えられる。そのため DMI 定数、磁界パルス幅による反転確率の変 化より DMI 定数が測定できることが考えられる。この手法による調査手順を以下に示す。まず 初期磁化を面直上向きで熱揺らぎによる平衡状態と設定する。これに、面直下向きの磁界パル スをかけ、磁化反転シミュレーションを反復して行い、反転確率を求める。さらに磁界の大き さを変えて、反転確率の磁界依存性を求め、得られた計算結果より DMI 定数を推定する。
第
6
章 反転磁界による
DMI
の測定法
第 6 章では反転磁界による DMI 定数の測定法について述べる。第 5 章で説明したように DMI により磁性体内部の磁化構造が変化するために、DMI 定数の測定に反転磁界を利用することが 可能であると考えられる。反転磁界は面直下向きと面内方向にそれぞれ外部磁界を印加するこ とで求めることができる面直反転磁界と面内反転磁界の 2 つがある。はじめに面直反転磁界と 面内反転磁界の計算方法について述べる。次に DMI 定数と円盤直径を変えて行った磁化反転シ ミュレーション結果について述べる。その後、得られた結果から DMI 定数による反転磁界への 影響の原因を調査する。最後に本論文における反転磁界による DMI 定数測定法を提案する。6.1
面直反転磁界の計算
まず面直 (z 軸) 方向に外部磁界をかけ、反転磁界を求めた。ここで求める反転磁界を「面直 反転磁界」と呼び、Hsw z と表す。DMI 定数と円盤直径による面直反転磁界の変化を調査した。6.1.1
計算条件
計算条件を以下に示す。 時間刻み : ∆t = 0.01 ps DMI定数 : D = 0.0, 0.2, 0.4, 0.6, 0.8, 1.0 erg/cm2 円盤直径 : d = 16, 20, 28, 36, 44, 52, 60, 68, 76, 84, 92, 100 nm 損失定数 : α = 1.0 収束条件 : tq(t) tq0 < 10 −4(tq :磁気モーメントのトルク) 反転条件 : hmzi < −0.956.1.2
計算結果
まず面直反転磁界を求めるために、外部磁界を印加し、面直磁化の時間変化をシミュレーショ ンにより求めた。磁化反転シミュレーションの例として、直径サイズ d = 52, 60 ns での結果を 図 6.1,6.2 に示す (D = 0.0, 1.0 erg/cm2)。-1 -0.5 0 0.5 1 0 0.5 1 1.5 2 2.5 3 <m z > t(ns) HzEXT=1100 Oe HzEXT=1700 Oe (a) D = 0.0 erg/cm2 -1 -0.5 0 0.5 1 0 0.5 1 1.5 2 2.5 3 <m z > t(ns) HzEXT=1100 Oe HzEXT=1700 Oe (b) D = 1.0 erg/cm2 図 6.1: d = 52 nm での外部磁界による面直磁化の時間変化 -1 -0.5 0 0.5 1 0 2 4 6 8 10 12 14 <m z > t(ns) HzEXT=1100 Oe HzEXT=1216 Oe (a) D = 0.0 erg/cm2 -1 -0.5 0 0.5 1 0 5 10 15 20 <m z > t(ns) HzEXT=1100 Oe HzEXT=1216 Oe (b) D = 1.0 erg/cm2 図 6.2: d = 60 nm での外部磁界による面直磁化の時間変化 図 6.1 では、D = 0.0, 1.0 erg/cm2共に HEXT z = 1100 Oeでは磁化反転せず、HzEXT = 1700 Oe では磁化反転することを表す。同様に図 6.2 より、D = 0.0, 1.0 erg/cm2共に HEXT z = 1100 Oe では磁化反転せず、HEXT z = 1216 Oeでは磁化反転することを表す。シミュレーションにより 求めた円盤直径 d 及び DMI 定数による面直反転磁界の変化を図 6.3 に示す。 面直反転磁界は、円盤直径の増加とともに減少するが、DMI による変化はあまりないことが わかった。円盤直径による面直反転磁界の減少の原因として、形状の変化による実効異方性磁 界の変化が考えられるため、後に議論を行う。
6.2
面内反転磁界の計算
次に、面内 (x 軸) 方向に外部磁界をかけ、反転磁界を求めた。ここで求める反転磁界を「面 内反転磁界」と呼び、Hsw x と表す。そこで DMI 定数と円盤直径による面直反転磁界の変化を調 査した。0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 Hz sw (kOe) d(nm) D=0.0 erg/cm2 0.2 erg/cm2 0.4 erg/cm2 0.6 erg/cm2 0.8 erg/cm2 1.0 erg/cm2 図 6.3: 円盤直径と DMI による面直反転磁界の変化
6.2.1
計算条件
計算条件を以下に示す。 時間刻み : ∆t = 0.01 ps DMI定数 : D = 0.0, 0.2, 0.4, 0.6, 0.8, 1.0 erg/cm2 円盤直径 : d = 16, 20, 28, 36, 44, 52, 60, 68, 76, 84, 92, 100 nm 損失定数 : α = 1.0 収束条件 : tq(t) tq0 < 10 −4(tq :磁気モーメントのトルク) 反転条件 : hmxi > 0.956.2.2
計算結果
まず面直反転磁界の計算と同様に面内反転磁界を求めるために、外部磁界を印加し、面内磁 化の時間変化をシミュレーションにより求めた。磁化反転シミュレーションの例として、直径 サイズ d = 60 ns での結果を図 6.4 に示す (D = 0.0, 1.0 erg/cm2)。図 6.4(a) では、D = 0.0 erg/cm2で HEXT
x = 900 Oeでは磁化反転せず、HzEXT = 1500 Oeで は磁化反転することを表す。これに対して、図 6.4(b) では、D = 0.0 erg/cm2で HEXT x = 900 Oe では磁化反転せず、、HEXT z = 1800 Oeでは磁化反転することを表す。シミュレーションにより 求めた円盤直径 d 及び DMI 定数による面内反転磁界の変化を図 6.5 に示す。 面内反転磁界は、面直反転磁界と同様に円盤直径の増加とともに変化するが、さらに DMI に より変化することがわかった。特に D = 1.0 erg/cm2では、d = 44 nm にて反転磁界が極小を 取ったのち、d = 68 nm まで反転磁界が増加した。このため DMI により、反転磁界は円盤直径 の増加により単調に減少するわけではないことがわかった。
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 <m x > t(ns) HxEXT= 900 Oe HxEXT=1500 Oe <mx>=0.95 (a) D = 0.0 erg/cm2 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 <m x > t(ns) HxEXT= 900 Oe HxEXT=1800 Oe <mx>=0.95 (b) D = 1.0 erg/cm2 図 6.4: d = 60 nm での外部磁界による面内磁化の時間変化 0 0.5 1 1.5 2 2.5 3 3.5 4 0 20 40 60 80 100 Hx sw (kOe) d(nm) D=0.0 erg/cm2 0.2 erg/cm2 0.4 erg/cm2 0.6 erg/cm2 0.8 erg/cm2 1.0 erg/cm2 図 6.5: 円盤直径と DMI による面内反転磁界の変化
6.3
面直、面内反転磁界の考察
6.3.1
実効異方性磁界との比較
6.1では面直反転磁界、6.2 では面内反転磁界の円盤直径および DMI 定数の依存性について 調査した。いずれの場合も、DMI 定数が 0.0 erg/cm2の場合、円盤直径の増加ととも反転磁界 が減少することがわかった。原因として円盤直径の変化による実効異方性磁界の変化が影響し ていると考えられる。静磁界を考慮した磁性体に実効的に働く実効異方性磁界は式 (6.1) で表 される [23]。 HKef f ≃ HK + HD (6.1) HKef f(Oe) :実効異方性磁界 図 6.6 に、円盤直径による「面直及び面内反転磁界」と「実効異方性磁界」を示す。 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 Hz (kOe) d(nm) HKeff Hzsw Hxsw 図 6.6: 実効異方性磁界と面直及び面内反転磁界の比較 膜厚一定の場合、円盤直径の増加により静磁界が増大する。静磁界は磁気モーメントの向き と反平行の方向に働くため、静磁界の増大により実効異方性は減少する。これより円盤直径の 増加による反転磁界の減少は実効異方性磁界により説明される。6.3.2
DMI
磁界による反転磁界への影響
6.1の結果より、面直反転磁界は DMI による変化があまりないことがわかった。ここでは DMI 磁界による実効磁界への影響を調査するために、外部磁界の元で磁化反転する原子磁気モーメ ントに作用する DMI 磁界の時間変化に注目する。-1 -0.5 0 0.5 1 0 1 2 3 4 5 6 40 60 80 100 120 140 160 180 200 220 240 <m z > <H DMI > t(ns) <mz> <HDMI> (a) HEXT = 1100 Oe -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 0 1 2 3 4 5 6 0 100 200 300 400 500 600 <m z > <H DMI > t(ns) <mz> <HDMI> (b) HEXT = 1300 Oe 図 6.7: D = 1.0 erg/cm2での平均磁化と DMI 磁界の時間変化 -200 -100 0 100 200 300 400 500 0 1 2 3 4 5 6 H DMI (Oe) t(ns) <HxDMI> <HyDMI> <HzDMI> 図 6.8: D = 1.0 erg/cm2での DMI 磁界の各成分での時間変化 図 6.7 に、円盤直径 60 nm、外部磁界 1100, 1300 Oe での磁化の DMI 磁界の時間変化を示す。 図中の赤線は平均磁化の磁間変化 (左軸に対応)、緑線は DMI 磁界の時間変化 (右軸に対応) を
表す。外部磁界で反転しない場合 (HEXT = 1100 Oe(図 6.7(a))) に現れる DMI 磁界は約 220 Oe
であった。一方、外部磁界で反転する場合 (HEXT = 1300 Oe(図 6.7(b))) に現れる DMI 磁界は、
磁化反転前は約 220 Oe であるが、磁化反転中に DMI 磁界が約 550 Oe と大きくなり、磁化反転 後は約 150 Oe となった。図 6.3.2 に、DMI 磁界の各成分の時間変化を示す。図より、基本的に は DMI 磁界は z 成分のみに働き、x, y 成分については磁化反転中に作用することがわかった。 面直反転磁界が DMI により変化しない理由として、DMI 磁界の面直成分が外部磁界の大きさ によらずに一定であることが考えられる。
6.3.3
外部磁界による面内磁化の変化
ここでは面内反転磁界の DMI 定数依存性について考察する。面内反転磁界の場合、初期磁化 を面直上向きの飽和状態としているため、初期磁界状態では x 方向の磁化成分は 0 となる。+x 軸方向に外部磁界を加えることにより、面内方向に磁化が飽和される。ここで、外部磁界の大 きさにより磁化が変化する様子を表すグラフである初磁化曲線を利用する。円盤直径が 40,60 nmでの初磁化曲線を図 6.9 に示す。図中の mx = 0.95は反転基準 hmxi > 0.95 を表す。 0 0.2 0.4 0.6 0.8 1 1.2 0 500 1000 1500 2000 2500 3000 <m x > Hx(Oe) D=0.0 erg/cm2 D=1.0 erg/cm2 mx=0.95 (a) 円盤直径 40 nm 0 0.2 0.4 0.6 0.8 1 1.2 0 500 1000 1500 2000 <m x > Hx(Oe) D=0.0 erg/cm2 D=1.0 erg/cm2 mx=0.95 (b) 円盤直径 60 nm 図 6.9: 面内磁化の初磁化曲線 図より、DMI により初磁化曲線の形状が変化することがわかった。D = 0.0 erg/cm2では、 磁界により面内磁化は単調に増加した。これに対して、D = 1.0 erg/cm2では、d = 40, 60 nm でそれぞれ 1700,1100 Oe にて面内磁化が急激に増加した後、飽和状態となった。このように、 外部磁界による面内磁化の変化量が DMI 定数により変化するため、DMI により面内反転磁界 が変化することがわかった。6.3.4
面直、面内反転磁界の差
6.1では面直反転磁界、6.2 では面内反転磁界の円盤直径および DMI の依存性について調査 した。円盤直径の影響は両方の反転磁界で確認できたが、DMI の影響は面内反転磁界のみで確 認できた。図 6.10 に、DMI 定数と円盤直径による面直、面内反転磁界の変化を示す。 図 6.10 において、面直反転磁界と面内反転磁界が交差する円盤直径が少なくとも 1 つ存在する ことがわかった。面内及び面直反転磁界への DMI と直径の影響が異なることから、この 2 つの磁 界の差と DMI 定数の関係を調べた。2 つの反転磁界の差を ∆Hswと表記し、∆Hsw = Hsw x −Hzsw とした。図 6.11 に、円盤直径による ∆Hswの変化の結果を示す。 図 6.11 で、∆Hsw = 0 Oeとなる円盤直径は、図 6.10 で面直反転磁界と面内反転磁界が交差 する円盤直径と同義である。そのため、図 6.12 に ∆Hsw = 0 Oeとなる円盤直径を表す。図より DMIの増加と共により ∆Hsw = 0 Oeとなる円盤直径が増加することがわかった。「∆Hsw = 0 Oe となる円盤直径」が DMI によって変化するため、この直径を求めることにより DMI を測定す ることが可能であることがわかった。0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (a) D = 0.0 erg/cm2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (b) D = 0.2 erg/cm2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (c) D = 0.4 erg/cm2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (d) D = 0.6 erg/cm2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (e) D = 0.8 erg/cm2 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 20 40 60 80 100 H sw (kOe) d(nm) Hzsw Hxsw (f) D = 1.0 erg/cm2 図 6.10: 円盤直径による面直、面内反転磁界の比較
-0.2 0 0.2 0.4 0.6 0.8 1 0 20 40 60 80 100 ∆ H sw (kOe) d(nm) D=0.0 erg/cm2 0.2 erg/cm2 0.4 erg/cm2 0.6 erg/cm2 0.8 erg/cm2 1.0 erg/cm2 図 6.11: 各 DMI 定数での円盤直径による ∆Hsw の変化 0 10 20 30 40 50 0 0.2 0.4 0.6 0.8 1 d (nm) D(erg/cm2) 図 6.12: DMI 定数による ∆Hsw = 0 Oeとなる 円盤直径
6.4
まとめ
本章では、反転磁界による DMI 定数の測定法を検討するために、磁化反転シミュレーション より面直、面内反転磁界の 2 つの DMI 定数及び円盤直径依存性を調査した。その結果、円盤直 径増大により面直、面内反転磁界が減少することがわかった。また面直反転磁界は DMI による 変化がなく、面内反転磁界は DMI による変化が現れることがわかった。円盤直径増大による反 転磁界の減少は、実効磁気異方性により説明できることがわかった。また DMI による面内反転 磁界の変化は、DMI による面内磁化の初磁化曲線の変化により説明できることがわかった。面 内、面直反転磁界の差 ∆Hswを計算し、差が 0 になる円盤直径を求めたところ、DMI の増加に より ∆Hsw = 0 Oeとなる円盤直径が単調増加することがわかった。これらの結果より、本章 にて提案する反転磁界による DMI 定数測定法を以下に示す。 1. 円盤状の直径の異なる試料を複数用意する。 2. 「円盤状の測定試料」の磁化を面直上向きに飽和させる。(初期値) 3. 測定試料に対して、面直方向と面内方向に磁界をかけ、それぞれの反転磁界を測定する。 4. 3で求めた 2 つの反転磁界の差 ∆Hswを計算する。 5. ∆Hsw = 0 Oeとなる円盤直径より、DMI 定数を求める。第
7
章 反転確率による
DMI
の測定法
第 7 章では反転確率による DMI 定数の測定法について述べる。現実世界では熱エネルギーに より原子磁気モーメントが振動している。熱揺らぎのランダム性のために原子磁気モーメント の磁化反転は確率的に起こる [22]。このため熱揺らぎを考慮した磁化反転では反転磁界ではな く、磁界による反転確率の変化を求めることになる。ここで DMI が作用する磁性体を考える。 DMIが作用する磁性体では、DMI により磁性体内部の磁化構造が変化するため、反転確率が変 化することが考えられる。本章では、熱揺らぎを加えた磁化反転シミュレータを用いてシミュ レーションを行う。まず熱揺らぎを加えた磁化反転シミュレータの妥当性をマクロスピンモデ ルにて確認した。次にマイクロマグネティックモデルにおいて熱揺らぎにより計算の安定性が 変化するため、計算の安定性を保証するための条件を調べる。その後、得られた安定条件を用 いて反転確率を計算する。反転確率のシミュレーションは DMI 定数と磁界パルス幅を変えて行 い、反転確率の外部磁界依存性を求める。最後にシミュレーションにより得られた結果を元に 本論文における反転確率による DMI 定数測定法を提案する。7.1
シミュレータの妥当性の調査
第 6 章までは、熱揺らぎ非考慮の磁化反転シミュレータによるシミュレーションを行った。第 7章では、熱揺らぎを考慮した磁化反転シミュレータによるシミュレーションを行う。ここで は、まず熱揺らぎを加えた磁化反転シミュレータの妥当性を調査する。シミュレータの妥当性 を調査するために解析式による比較が可能であるマクロスピンモデルを用いる。検証方法は、 初期磁化を ~m = (0, 0, 1)として、熱揺らぎによる磁化の平衡状態を解析式と比較することであ る。前述の通り、熱揺らぎを加えた単磁区微粒子での平衡状態の磁化角度の解析式は式 (7.1) で 表される [23]。 θ = r 1 ∆ (7.1)7.1.1
計算条件
シミュレータの妥当性に調査において使用した計算条件を以下に示す。 時間刻み : ∆t = 0.01, 0.1, 1 ps 損失定数 : α = 1.0 初期磁化 : ~m = (0.0, 0.0, 1.0) 体積 : v = 1.78 × 10−19 cm3 温度 : T = 300 K(室温)7.1.2
計算結果
計算で得られた磁化角度の時間平均の変化を図 7.1 に示す。 図中の赤線は、磁化角度の解析値を表す。時間経過とともに平衡状態が磁化角度の解析値に 近付くことが確認できた。特に α ≥ 0.1 にて、t ≤ 10 ns で磁化角度が解析値に近づくことがわ かった (図 7.1(a),7.1(b))。しかし、α ≤ 0.001 では t ≤ 100 ns でも磁化角度が安定しなかった (図 7.1(d))。このように損失定数が小さい場合、平衡状態になりにくくなる。しかしながら時 間刻みが小さい場合は、平衡状態に到達しやすい結果となった。また時間刻みを変えた場合、 ∆t = 0.01 psの場合では最も計算が安定しやすく、短時間で磁化角度の解析値に近付く結果と なった。これに対して ∆t = 1 ps の場合では、磁化角度の収束までに時間を最も多く要する結 果となった。シミュレーションにより求めた結果を元に、表 7.1.2 に熱揺らぎを考慮したマクロ スピンモデルでの損失定数と時間刻みによる計算の安定性を表す。 表 7.1: 熱揺らぎを考慮したマクロスピンモデルでの損失定数と時間刻みによる計算の安定性 α\∆t 1 ps 0.1 ps 0.01 ps 1.0 ○ ○ ○ 0.1 ○ ○ ○ 0.01 ○ ○ ○ 0.001 × △ △ 0.0001 × × ×7.2
計算の安定性の調査
7.1では、熱揺らぎを加えた磁化反転シミュレータの妥当性をマクロスピンモデルにて調査し た。ここからはマイクロマグネティックモデルにて熱揺らぎを加えた磁化反転シミュレータを 用いる。磁化反転シミュレータに熱揺らぎの項を加えた場合、時間刻みと損失定数の組合せに より計算が不安定になる場合があるため、計算が安定する条件を見つける必要がある。ここで は、使用する損失定数に合った時間刻みを調べた。7.1 の結果では、時間刻みが小さい場合に平 衡状態に到達しやすい結果となった。しかし、時間刻みを小さくしすぎた場合、計算時間を多 く要する。そのため、計算が安定する条件での最大の時間刻みを求める必要がある。7.2.1
計算条件
計算の安定性の調査に用いた計算条件を以下に示す。 時間刻み : ∆t = 0.02, 0.025, 0.05 ps 損失定数 : α = 1.0, 0.1 温度 : T = 300 K(室温)0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0 2 4 6 8 10 < θ > t(ns) sqrt(1/∆) ∆t=1 ps ∆t=0.1 ps ∆t=0.01 ps (a) α = 1.0 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0 2 4 6 8 10 < θ > t(ns) sqrt(1/∆) ∆t=1 ps ∆t=0.1 ps ∆t=0.01 ps (b) α = 0.1 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0 2 4 6 8 10 < θ > t(ns) sqrt(1/∆) ∆t=1 ps ∆t=0.1 ps ∆t=0.01 ps (c) α = 0.01 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0 10 20 30 40 50 60 70 80 90 100 < θ > t(ns) sqrt(1/∆) ∆t=1 ps ∆t=0.1 ps ∆t=0.01 ps (d) α = 0.001 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0 10 20 30 40 50 60 70 80 90 100 < θ > t(ns) sqrt(1/∆) ∆t=1 ps ∆t=0.1 ps ∆t=0.01 ps (e) α = 0.0001 図 7.1: マクロスピンモデルでの磁化角度の時間変化
0.88 0.9 0.92 0.94 0.96 0.98 1 0 2 4 6 8 10 <m z > t(ns) ∆t=0.02 ps ∆t=0.025 ps ∆t=0.05 ps ∆t=0.1 ps ∆t=0.5 ps (a) α = 1.0 0.92 0.93 0.94 0.95 0.96 0.97 0.98 0 2 4 6 8 10 <m z > t(ns) ∆t=0.02 ps ∆t=0.025 ps ∆t=0.05 ps ∆t=0.1 ps (b) α = 0.1 図 7.2: D = 0.0 erg/cm2平均磁化の時間変化 0.92 0.925 0.93 0.935 0.94 0.945 0.95 0.955 0.96 0.965 0.97 0 2 4 6 8 10 <m z > t(ns) ∆t=0.02 ps ∆t=0.025 ps ∆t=0.05 ps ∆t=0.1 ps (a) α = 1.0 0.92 0.93 0.94 0.95 0.96 0.97 0.98 0.99 0 2 4 6 8 10 <m z > t(ns) ∆t=0.02 ps ∆t=0.025 ps ∆t=0.05 ps ∆t=0.1 ps (b) α = 0.1 図 7.3: D = 1.0 erg/cm2平均磁化の時間変化
![図 2.2: 磁気モーメントの概念図 2.3 LLG 方程式 LLG(Landau-Lifshitz-Gilbert) 方程式とは、磁界による原子磁気モーメントの運動を表す方 程式である。LLG 方程式は、式 (2.1) で表される [8]。 ~˙ M = −| γ | ( M~ × H) +~ α M ( M~ × M~˙ ) (2.1) γ : 磁気回転比 α : 損失定数 M~ : 原子磁気モーメント H~ : 磁界 式 (2.1) において第 1 項をジャイロ項と呼び、磁界周りの磁気モーメントの歳差](https://thumb-ap.123doks.com/thumbv2/123deta/7725203.1711294/9.918.343.551.65.287/モーメントモーメントモーメントにおいてジャイロモーメント.webp)
