L. W. Conner and A. P. Malozemoff
physical review B 43(1)402,1991
概略
我々はディスク型の超伝導体の臨界状態における磁場分布と磁化の値を数値計算によって示す。臨界電流の磁場依存性にKimモデルを使って、以下のことを示す。これまでのflux line の湾曲を無視するモデルではゼロ磁場でのremanent moment が半径に対して(on radius)非線形に依存するが、湾曲を考慮した場合数値計算によるとほとんど直線の依存性となって、Beanモデルと比べると小さめの傾きになる。近似的にこの減少した傾きを表す式を求めた。
イントロ
この論文では臨界状態におけるディスク型の超伝導体の磁場分布と磁気モーメントの数値計算を示す。......巨大磁化効果によってこの状態(this configuration)での解釈が複雑になっているが、このケースを理解することは重要である。なぜならalignment problems and sigunal strength problemsがin-plain測定を難しくしているし、それ故異方的な物質ではいつも問題になる。
通常の測定では形状の効果（有限の大きさ）を取り扱うために、反磁性係数Ｎを用いるがこれはサンプルの形状を楕円体近似を使うもので、また加えて印加磁場に-N*Mだけ修正をする。ここでMは磁化である。......しかしこれは定性的にも定量的にも誤りであると指摘されている。DaeumlingとLarbalestier は計算によってＲ（半径）がＤ（厚さ）に比べて大きいディスクサンプルで均一に流れる同心円上のJcが湾曲した磁束線を与えることを示したが、そのときでも通常のBeanモデルの、残留磁化Mremは
Mrem=Jc*R/3
で、変化しない。しかし臨界電流密度JcはKimモデルに従う。DaeumlingとLarbalestier は磁化が磁場と半径に複雑な形で依存していることを示した。彼等は磁場に依存した磁化Mremを一例計算している。
ここで私達は残留磁化の半径依存性に注目する。式１は広く臨界電流の決定に用いられているから、半径依存性を実験的に求めることで式１の正確さを検討することは重要だ。現在の所実験結果は曖昧だ。....
線形の結果はBeanモデルとバルク電流を(Edge currentではなく）支持していると解釈される。しかし通常の直線状の磁束量子線を考えると矛盾する。なぜならこの場合は磁化Ｍは２章に示すように、Jcが磁場に依存しようと非線形だからである。この矛盾を解決するにはvortexに湾曲を取り入れればよい。......我々はこの計算を３章で示す。最後の４章で結果を解釈して、この複雑な振る舞いに物理的洞察を与える近似式を示すが、それによってdMrem/dRがBeanモデルよりも小さくなることが示せる。ここではJcの磁場依存性を等方的として考えている。異方性を考慮するのは後の仕事でやりたい。
2.conventional model
我々は通常のモデルでJcに磁場依存性を取り入れたいわゆるKimモデルを用いて磁気モーメントの半径依存性を計算した。通常のモデルとは面に垂直に磁束線がなっている場合である。マックスウェル方程式から
rot(H)=J (3)
一般的にはdHz/drとdHr/dzがゼロではないが磁束量子が面に垂直であるとすれば、dHz/dr=Hr=0で(3)式はdHz/dr=Jになる。.....中略...
Jcの磁場依存性を考慮して(3)式を積分すると下記、
(4)式
次に、磁気モーメントMを無視できるぐらい小さいとして、Bをu0Hに等しいとみなす。これはHc1を他の磁場に比べて小さいと考えるのに等しい（たとえば式(2)の中のH0、あるいは自己磁場と比較して）。Bを体積全体で積分して体積で割れば（体積あたりにすると）磁化が求まる。
結果を示すために以下の規格化を行う（この論文では以下この変数を用いる）。....中略....
(5),(6),(7),(8)式
この表記法を用いると規格化されたmは以下の式になる。
(9)式
これをFig.1に実線で示す。明らかに非線形である。この非線形性はH(r)からきているがこれは端から中心にかけて大きくなってJcの減少をもたらしている。明らかにRが小さければこのセルフフィールドはH0に比較して小さいから、Jcを定数としてよい。極限において式(9)はBeanの式に一致する。これはm=h0*x/3であるが、これをドットでFig.1に示す。
YBCOの薄膜と単結晶を用いた典型的な実験はいかに示すように非線形領域になる。低温でJcは単結晶で2x10^10[A/m^2]で（薄膜はそれ以上）式(2)にfitさせればu0*Hは0.5Teslaでu0*alphaは10^10[AT/m^2]である。自己磁場がH0に匹敵する半径はH0/Jc、あるいは25umとなって、ほとんどのサンプルの寸法よりも小さい。
それなら残留磁化の半径に対する線形性はなぜだろう。答えは問題を完全な形で取り扱えば得られるが、特に磁場は垂直成分だけでなく湾曲していて、dHz/drよりもdHr/dzを作り出していることによる。
3.Numerical Calculations
我々はここでバルク電流の流れる超伝導体における磁場と電流分布の計算を記述する。ここでは侵入長よりも厚さが十分に大きいと仮定してそれ故にfluxの湾曲やLondon表面電流によるエネルギーを無視する。円筒状である場合半径aの位置にあるリング状の電流が作る磁場は式(10),(11)...中略...それぞれのpointの磁場はこれを積分して得られる。Hc1が小さいと仮定したからこのサンプル内における磁場をvortexによって作られた空間的に平均されたlocal磁場と考える。（意味が不明）
式(10)(11)を積分するのにRomberg法と似た方法を採ったが拡張された台形法(trapezpidrule)を使う代わりに拡張されたミッドポイント法を用いた。これには２つの理由がある。数値積分とは違って両方の方法ともステップサイズの偶数次の項によるエラー級数を有す（本当か？）。これにより積分区間を細かくしていけば急速に誤差は小さくなるが台形法とは違ってミッドポイント法はsubintervalにおいて積分の評価を必要としない。これによりa=rにおける積分の発散を避けることができる。Romberg法によれば積分値の最初の見直し(refinements)はステップサイズを小さくすることでなされるが積分値を再評価して前の値と組み合わせてエラーを減少させる。これはrとzの２重積分であるからrに関する積分はそれぞれのzの値のsubintervalにおいて必要な精度で計算されなければならない。それ故にzのステップサイズが小さくなって精度が上がるにつれてrに関しても増加しなければならない。計算にかかる時間を減らすために、外挿のテクニックを合わせて用いてステップサイズを小さくする。それぞれのステップサイズの積分値を記録する。この値を用いてnevil's algorithm の内挿多項式を作りゼロステップの時（極限値）の積分値を外挿する。もし誤差の見積もりが必要な精度に達しなければステップサイズは小さくなり積分値は拡張されたミッドポイント法により再評価され、内挿多項式は新しく計算し直されて再びゼロステップでの値が求められる。
この方法を用いてディスク内のHrとHzを計算した。電流値が一定の場合Daeumling and Larbalestierの結果と同じものを得た。それから等方的ではあるが磁場に依存した電流密度（臨界電流密度）の場合(Kim model）を計算した。磁場の強さは電流密度によりまた電流密度は磁場の強さに依存するのでこの等式は繰り返し解かれる。（フィードバックが必要だ）あるポジションの磁場はディスク全体の電流密度に依存するので各地点の磁場を同時に求めておくことが必要である。そうするために計算は電流密度が一定かあるいは実験から求めた（磁場依存性の）値を用いて始める。グリッドを決定し電流の初期値を定めて式(10)(11)をそれぞれのグリッドで積分して磁場を求める。ディスク全体にわたって電流密度を式(2)を用いて計算する。その電流から新しい磁場の値を計算する。この手順を磁場と電流の値がそれぞれのグリッドで矛盾しなくなるまで続ける。HrとHzの積分をする時にはグリッド上ではないところの電流値を決定することが必要である。２次元的に内挿する事が必要であるから潜在的に意味のある（無視できない）誤差が生じて電流密度は時として急激な変化を見せるが、我々は内挿多項式を用いて電流値を計算した。再びNevillesのアルゴリズムを用いてあるポイントrの電流値を内挿するために３から５ポイントの値を用いた。（不明なところが？）それから必要なポイントの電流密度を求めるのにこれらの値から内挿した。
収束にはいくつかのファクターが影響する。通常の問題である発散、振動、カオス的振る舞いなどに加えてこの等式固有の問題が現れる。HrとHzはJによっているがJもHrとHzによっている。この相互依存性HrとHrのどちらにエラーが現れても次の回の計算に影響する。Hrの計算においてどのポイントでもすぐに収束するかもしれない。しかしHrの誤差は減少するにもかかわらずHzの誤差は増加するかもしれない。結局Hzの誤差によってHrが正確な値からずれてしまうだろう。計算の途中をよくモニターして実際に解に近づいているかどうか確認する事が必要である。他に重要な問題として誤差の伝搬がある。これはHrもHzもディスク全体にわたる積分であることによる。それ故あるポイントにおける磁場の値の誤差が他のポイントの磁場の値に影響する。これによって積分に必要とされる精度は収束するかどうかに大きく影響してくる。....中略.....我々はちょうどよい妥協点として積分における誤差を10^-3にした。収束に影響する他の要素としては磁場を計算するポイントの数がある。計算するポイントに数が増えれば明らかにそれぞれの回での計算に要する時間が増加するが、収束に要する繰り返し回数は減少する。これにはいくつかの理由が考えられるが、一つは電流密度がそれほど場所によって変化しないことがあげられる。ポイント間で内挿をなめらかに変化させれば（意味が不明）式(10)と式(11)を性格に積分する事が可能になって、それにより今度はHrとHzの積分値が正確になる。我々は能率が最適化されたポイント数を見つけることはできなっかた。新しいパラメータを用いるときには推測や経験から、計算の進行ぐあいをみてポイント数を減らしたり増やしたりした。普通は2000から7000ポイントを用いた。収束を早めるためにいくつかのテクニックを用いたがその中に、前回の計算値を用いた重み付き平均、選択的平均（不明）、何ポイントかの点を用いたデータのスムージングがある。
収束に関して最後に、それが唯一の数学的解ではないことを指摘したい。我々がこの計算結果に求めるもののすべては得られた電流値が式(2)(10)(11)を通じてそれ自身を矛盾なく作り出すことである。物理的に唯一の解があると信じる一方で数学的要求だけを満たすいろいろな電流分布があるかもしれない。妥当性を高めるには（信頼性をあげるには）異なる電流分布やポイント数の増加などを試みてみて正しい解を得る必要がある。（どんな方法でやっても同じ結果を得る。）
例をFig2aにあげる。ディスク表面におけるHrとディスクのmid plainでのHzを描いている。R/Dの計算にはHzにいくらかばらつきが見られる。R/Dが大きいほど計算においてたくさんのポイントをとる必要がある。この場合には6480ポイントを用いた。もっと小さな半径での計算結果によればポイントを増やせばもっとスムーズになるが特に重要な変化はなかった。十分に計算には時間がかかっていたので、これ以上重要な情報がでてきそうにはなかったのでR/Dが15の場合はこれ以上ポイントを増やさなかった。
結果はDaeumlingとLarbalestierのものやFrankelのものと定性的には同じである。Fig.2において特徴的なことは表面でのHrの大きさである。対称性からいえばHrは中心でゼロであり下面でその符号を変える。Fig.2(a)によれば（ｂの間違いでは？）Hrの値はｚ依存性が強いが（ディスク厚方向）半径方向に関しては中心とエッジの付近をのぞいて依存性が弱く、半径に対しても強くは依存しない。これはfig.2(b)で確認できるが、この図はHrとHzの厚さ方向の依存性を示したものである。
fig.2(b)はまたHzがｚ方向でほぼ一定であることを示している。（a）図によればｒ方向には分布がある。この傾きは図の真ん中あたりで接線をひいてみると見積もれる。D/Rの比による垂直方向の傾きからは（たとえばJc一定の場合の値からは）小さな値となっている。（意味が不明）
Jcの分布を等高線で表した図からはもっと多くのことが得られる。(fig.3)これはJcはKimの形式を用いて、h0は0.5、R/D=5である。等高線は0.6（サンプルの中心）から2（サンプルの端付近）まで0.05刻みで書いた。ピークがこのあたりにある理由はfig.2から分かるように、ここでHz、Hrがゼロであるからでそれ故Jcが大きい。（本当か？）
上記の結論と同じ内容のことがすでに文献５，６で述べられている。それによれば磁場の傾きを決める主因は一般的に考えられているようにdHz/drではなくてdHr/dzであり、これはR/Dに比例して膜や単結晶プレートでは非常に大きい。これはvortexが強く湾曲していることを示し、臨界状態の傾きが厚さ方向で支配的であることを示す。（？）Jcが一定の場合との主要な違いは磁場の傾きの絶対値である。dHr/dzがJcに等しくなるというより、低磁場領域での値よりも小さな値をとるが、これは明らかに自己磁場がH0程度の大きさとなってJcを小さくしているからである。
最後に試料全体での体積あたりの磁気モーメントを求める。計算は基格化半径が２から200まで、基格化磁場が0.1から2までの範囲で行った。
結果をfig.1に示す(H0=1)。これを通常のモデル、dHz/drモデルと比較した。計算値は丸で表している。線は以前のモデル？より上に位置していて半径に対してかなり線形依存である。明らかに曲線はオリジンの正確な外挿ではない。なぜなら半径が小さい極限では薄いシリンダー状の形状と見なせるので、通常のモデルの方が正確なものとなり、傾きは式１に近づき、そのときはx/3となる。それにも関わらずこの訂正は小さい。直線はゼロでない２点のデータから（異なるh0から計算によって求めた値から）決定されるが、それをfig.4に示す。
4.discussion
計算によって求めた磁化が近似的に、半径に対して線形依存であることがこの論文の主な結論である。磁場に依存する電流密度を仮定した通常のモデルが半径に対して非線形依存をするのとは対称的である。この元になる理由として臨界状態の傾きが半径よりむしろ厚さで決まっていることである。この見解は半径に対する線形依存が磁場に強く依存するJcの場合に見受けられることを説明するかもしれない。
この線形性を決定しているのは何かという疑問が残るが、それは単純なBeanモデルが予言するものでもある。（ここではJc,0は低磁場極限をとった。）これからすると詳細な計算をすることなく残留磁場の測定からJcを導くことができそうである。先のDaeumlingとLarbalestierの結果とFig.2の類似性は我々の結果に単純な説明が可能であることと、解析的公式化を導き出すときの問題点を示している。キーとなる点は磁場の傾きが主に厚さ方向を向いていることである。（fig.3の等高線が大部分垂直であることを指すと思う）また厚さは薄いので厚さに関して電流を平均することが（厚さ方向で一定の電流値とみなす）できて、この平均された電流値はほとんどの半径において同じであり（fig.3で等高線の間隔がだいたい一定であると見なせることを指す）主な性質を決定していると十分見なせる。fig.3は周辺付近をのぞいては十分に均一であり、この考えが近似の最初として妥当であることを示している。
この計算をするために我々は式３のdHz/drを無視した。（これは２章の反対の極限である）これはD/Rの程度で正しい。こうしてHrの厚さ依存性を計算することができた。得られた式は式４と一致して、ただしzをR-rで置き換える。この結果をKimの式に代入して厚さ方向で積分すると平均された電流値に関して式１２を得る。
ここでh0=H0/(alpha*D)^(1/2)である。この平均された電流は半径に依存せず均一であると考えられ、大きさは式１であたえられるが、ただしこのときのJcはJc,avである。Jc,0=alpha/H0を使って傾き減少係数なるものを定義できる。
式13
fig.1によればこの傾き減少係数は本当の傾きとBeanモデルとの比を表す。式13をh0の関数としてfig.4に書いたが、以前の計算から求めた値と比較できる。その一致は十分なもので、物理的概念の正確さを示している。式１３はh0が１以上での飽和を示しているが、これはH0=Jc,0Dの条件に対応する。その意味はH0が自己磁場のJc,0D程度になれば式３の電流密度は磁場に依らなくなり、そのため傾き減少係数は１に近づく。その考えは(The formula)fig.1においてxが１以上の領域で正しい（xは式８の規格化半径）これはJc,0R>H0に等しい条件である。
まとめると、式１３はBeanモデルと残留磁化からの臨界電流密度の決定に関する単純な補正を与える。これから先の研究では、残留磁化の外部磁場依存性や臨界電流の異方性などが興味深い。いずれもここでは無視している。これらは今後の研究課題である。