令和５年度

修士論文

超伝導リザバーの精度向上に関する研究

九州工業大学情報工学府
情報創成工学専攻 物理情報工学専門分野

学生番号 226E0301

有田 拳
指導教員：小田部 荘司

目次
目次 .................................................................................................................................. 0
第1章

序論 ................................................................................................................ 1

1.1

はじめに ......................................................................................................................1

1.2

超伝導現象 ..................................................................................................................2

1.2.1

超伝導現象の概要 ..............................................................................................2

1.2.2

超伝導体の分類 ..................................................................................................2

1.3

Time-Dependent Ginzburg-Landau 方程式 .................................................................3

1.3.1

Ginzburg-Landau 方程式 .....................................................................................3

1.3.2

Time-Dependent Ginzburg-Landau 方程式..........................................................4

1.4磁束の量子化 ..............................................................................................................5

1.5

磁束ピンニング機構 ..................................................................................................5

1.5.1

要素的ピン力 ......................................................................................................6

1.5.2

巨視的ピン力密度 ..............................................................................................6

1.6

Affine Integrator(AFI) ................................................................................................7

1.6.1

1 次元系における空間・時間離散化................................................................8

1.6.2

2 次元系での AFI .............................................................................................. 11

1.7

リザバーコンピューティング ................................................................................14

1.7.1

リザバーコンピューティングの概念 .............................................................14

1.7.2

エコーステートネットワーク.........................................................................15

1.7.3

バッチ学習(線形回帰とリッジ回帰) ..............................................................16

1.7.4

評価指標 ............................................................................................................18

1.7.5

NARMA モデル ................................................................................................181.7.6

非線形-メモリタスク .......................................................................................19

1.8

研究目的 ....................................................................................................................19

第2章

実装および計算方法 .................................................................................. 20

2.1

実装する TDGL 方程式と固定パラメータの規格化 ...............................................20

2.2

空間離散化 ................................................................................................................212.3

境界条件の処理 ........................................................................................................23

2.3.1

リンク変数の境界条件 ....................................................................................24

2.3.2

オーダーパラメータの境界条件.....................................................................26

2.4

初期条件と時間離散化 ............................................................................................27

2.5

電界、磁界、電流密度の実装 ................................................................................28

2.5.1

電界 ....................................................................................................................28

2.5.2

磁界 ....................................................................................................................292.5.3

電流密度 ............................................................................................................292.6

ピンの導入 ................................................................................................................30

2.7

描画............................................................................................................................30

2.8

リザバーコンピューティング ................................................................................30

2.8.1

ピンの数の変更 ................................................................................................322.8.2ピンの配置の変更 ............................................................................................32

2.8.3

ピン 1 つあたりのピン力の変更.....................................................................32

第3章

結果と考察 .................................................................................................. 33

3.1

量子化磁束および電磁現象の可視化.....................................................................33

3.1.1

ピンの数を変化させた場合 ............................................................................333.1.2ピンの配置を変化させた場合.........................................................................35

3.1.3

ピン 1 つあたりのピン力を変化させた場合 .................................................37

3.2

NARMA2 タスク ......................................................................................................39

3.2.1

ピンの数を変化させた場合 ............................................................................393.2.2ピンの配置を変化させた場合.........................................................................41

3.2.3

ピン 1 つあたりのピン力を変化させた場合 .................................................423.3

非線形-メモリタスク ...............................................................................................443.3.1

ピンの数を変化させた場合 ............................................................................44

3.3.2ピンの配置を変化させた場合.........................................................................47

3.3.3

ピン 1 つあたりのピン力を変化させた場合 .................................................50

第4章

まとめ .......................................................................................................... 53

参考文献 ........................................................................................................................ 54
謝辞 ................................................................................................................................ 57
研究業績 ........................................................................................................................ 58

図目次
Fig. 1.1 外部磁界に対する超伝導体の内部の磁束密度 .................................................3
Fig. 1.2 複数の磁束線のピン止め[25] .............................................................................7
Fig. 1.3 格子・リンク変数[29] .........................................................................................8
Fig. 1.4 2 次元系におけるダブルインデックスの定義(周期的境界条件)[29] ............12
リザバーコンピューティングの概念図を Fig. 1.5 に示す。英語の "reservoir "は液体を入れる容器、または日常使用する水を貯める貯水池を意味する。例えば、水面（リザバー）に石を落とすと波紋が広がる。この波紋は、投入された石の大きさや形によって変化する。つまり、水面の波紋は入力された情報を反映していると考えられる。 ................................................................................................14
Fig. 1.6 リザバーコンピューティングの概念図 ..........................................................14
Fig. 1.7 エコーステートネットワークの基本モデル...................................................15
Fig. 2.1 演算∇ ×のリンク変数による差分化およびループ変数の定義[29] ..............21
Fig. 2.2 境界付近の格子点およびリンク変数の配置[29] ............................................23
Fig. 2.3 東西境界における印加磁界、印加電流のリンク変数への反映[29].............25
Fig. 2.4 境界条件と基礎方程式によるオーダーパラメータの更新[29] ....................26
Fig. 2.5 電流の印加方法と電界の抽出方法の概略図...................................................31
Fig. 3.1 ピンの数ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。上から順に、ピンの数が 0、10、20、30、40、50 となっている。 .....................35
Fig. 3.2 ピンの配置ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。
上から順に、ピンを中央付近に配置した場合、左右境界付近に配置した場合、
上下境界付近に配置した場合を表している。 .....................................................37
Fig. 3.3 ピン 1 つあたりのピン力に関係するパラメータ𝛼ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。上から順に、𝛼 = 0, 5, 10, 15, 20のときを表している。 ............................................................................................................39
Fig. 3.4 各ピンの数における NARMA2 タスクの結果。横軸が時間、縦軸が出力を表している。(a)はピンの数が 0、(b)はピンの数が 10、(c)はピンの数が 20、(d)はピンの数が 30、(e)はピンの数が 40、(f)はピンの数が 50 のときの結果を表している。 ....................................................................................................................40
Fig. 3.5 電流密度の乱数波入力に対する電界の出力応答 ...........................................41
Fig. 3.6 各ピンの配置における NARMA2 タスクの結果。(a)がピンを中央付近に配置した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付近に配置した
場合の結果である。 ................................................................................................41
Fig. 3.7 各𝛼の値に対する NARMA2 タスクの結果。(a)が𝛼 = 0、(b)が𝛼 = 5、(c)が
𝛼 = 10、(d)が𝛼 = 15、(e)が𝛼 = 20に対応している。........................................43
Fig. 3.8 入出力データの FFT ..........................................................................................44
Fig. 3.9 各ピンの数における非線形-メモリタスクの R2 の𝜈依存性。(a)がピン 0 個、(b)がピン 10 個、(c)がピン 20 個、(d)がピン 30 個、(e)がピン 40 個、(f)がピン
50 個のときの結果を表している。........................................................................45
Fig. 3.10 各ピンの数における非線形-メモリタスクの R2 の𝜏依存性。(a)がピン 0 個、(b)がピン 10 個、(c)がピン 20 個、(d)がピン 30 個、(e)がピン 40 個、(f)がピン
50 個のときの結果を表している。........................................................................47
Fig. 3.11 各ピンの配置における非線形-メモリタスクの R2 の𝜈依存性。(a)がピンを中央付近に配置した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付
近に配置した場合の結果を表している。 .............................................................48
Fig. 3.12 各ピンの配置における非線形-メモリタスクの R2 の𝜏依存性。(a)がピンを中央付近に配置した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付
近に配置した場合の結果を表している。 .............................................................49
Fig. 3.13 各𝛼の値における非線形-メモリタスクの R2 の𝜈依存性。(a)が𝛼 = 0の場合、(b)が𝛼 = 5の場合、(c)が𝛼 = 10の場合、(d)が𝛼 = 15の場合、(e)が𝛼 = 20の場合の結果を表している。 ............................................................................................51
Fig. 3.14 各𝛼の値における非線形-メモリタスクの R2 の𝜏依存性。(a)が𝛼 = 0の場合、(b)が𝛼 = 5の場合、(c)が𝛼 = 10の場合、(d)が𝛼 = 15の場合、(e)が𝛼 = 20の場合の結果を表している。 ............................................................................................52

表目次
Table. 3.1 ピンの配置に対する R2..................................................................................42

第1章

序論

1.1 はじめに
近年、人工知能(AI)、中でも機械学習に関する研究が盛んに行われている[1-3]。その主流となるのがニューラルネットワークである[4]。従来のニューラルネットワークには、学習に膨大な時間とコストがかかるという欠点がある。リザバーコンピューティン
グの登場により学習にかかる時間とコストは大幅に削減された。リザバーコンピューテ
ィングとは、時系列データ処理に適した計算フレームワークである[5]。これはエコース
テートネットワーク[6-8]やリキッドステートマシン[9, 10]などのリカレントニューラル
ネットワークから派生したものである。リザバーコンピューティングは、入力を高次元空間にマッピングし、その高次元状態からパターン分析のための読み出しを行うリザバーで構成されている。リザバーの特性は固定されており、読み出しだけが線形回帰や分類などの単純な方法で学習される。したがって、他のリカレントニューラルネットワー
クと比較して、リザバーコンピューティングの主な利点は、迅速な学習と低い学習コストである。さらに、適応的な更新を行わないリザバーは、様々な物理システム、ボード、デバイスを使用したハードウェア実装に適している[11]。実際、物理リザバーコンピュ
ーティングは、スピントロニクス、神経細胞培養、ナノ磁性デバイスなど、多様な研究分野で注目を集めている[12-17]。最近では、量子力学に基づくリザバーコンピューティ
ングの研究も行われている[18-20]。本研究は、数ある非線形性を持つ物理現象の中でも、第 2 種超伝導体における量子化された磁束線の動きに関連する電流-電圧特性に注目して行われた。超伝導とは、極低温下において完全反磁性、電気抵抗ゼロの 2 つの特徴をもつ物理現象である。中でもピンニングセンターを含む第Ⅱ種超伝導体は、強い非線形性と記憶特性を持つため、リザバーに適していると考えられる。この超伝導現象を利用した物理リザバーを以下超伝導リザバーと呼ぶこととする。超伝導リザバーの有用性については、過去卒業論文にていくつかのリザバーコンピューティングのタスクを実行するという形で示した[21]。しかし、実際にリザバーとして運用できるほど高い精度を示すことはできなかった。これに関連するリザバーコンピューティングそのものの問題として、リザバー層はいわゆるブラックボックス化されており、リザバーの中身と時系列予測の結果の関係が読み取りづらいというものが挙げられる。リザバーコンピューティ
ングにおいて内部状態と予測精度の関係を明らかにすることは重要な課題である。

- 1 -

1.2

超伝導現象

超伝導とは、1911 年に Kamerlingh-Onnes により発見された現象である。液体ヘリウ
ムを用いた極低温における水銀の電気抵抗の測定中、絶対温度 4.2 K で急に測定不能なほど電気抵抗が小さくなったことから超伝導が発見された[22]。

1.2.1 超伝導現象の概要
ある温度以下で電気抵抗がゼロになる現象を超伝導といい、この状態を超伝導状態という。逆に電気抵抗がある状態を常伝導状態という。電気抵抗がゼロになる温度を臨界
温度𝑇c という。一般に、臨界温度は金属系超伝導体では低く、酸化物系超伝導体では高い。
超伝導のもう 1 つの特徴は、完全反磁性であることである。完全反磁性とは、磁場を印加しても内部に磁場が侵入しないこと、あるいは常伝導状態で磁場を印加した後、物質が冷えて超伝導状態になると内部の磁場が外部にはじき出されることを意味する。この特徴はマイスナー効果としても知られている。
超伝導状態になる条件は、𝑇c のほかに、臨界磁場𝐻c と臨界電流密度𝐽c がある。温度、
外部磁場、内部を流れる電流密度がそれぞれ𝑇c 、𝐻c 、𝐽c 以下のときに超伝導状態が実現する。

1.2.2 超伝導体の分類
超電導体は 2 つのタイプに分類される。これらは第 1 種超伝導体と第 2 種超伝導体と呼ばれている。これらは磁気特性によって区別される。一般に、磁束密度𝐵は、
𝐵 = 𝜇0 𝐻e + 𝜇0 𝑀

(1.1)

と表される。ここで𝐻e は外部磁界、𝑀は磁化、𝜇0 は真空の透磁率である。第 1 種超伝導体の超伝導体内部の磁束密度は Fig. 1.1(a)のように表される。𝐻e < 𝐻c では、𝑀 = −𝐻e となり、外部から磁場が与えられても内部の磁束密度はゼロのままである。この状態を超伝導状態と呼ぶ。𝐻e > 𝐻c では、𝑀 = 0となり、内部の磁束密度は、
𝐵 = 𝜇0 𝐻e(1.2)となる。このとき超伝導状態は壊れ常伝導状態へと相転移する。
第 2 種超伝導体の内部の磁束密度は Fig. 1.1(b)のように表される。第 2 種超伝導体では、下部臨界磁界𝐻c1 と上部臨界磁界𝐻c2 が存在し、𝐻e < 𝐻c1 では第 1 種超伝導体と同様に、内部の磁束密度はゼロのままである。𝐻c1 < 𝐻e < 𝐻c2 のとき、磁束は量子化された状態で超伝導体の内部に入る。そのため内部の磁束密度は徐々に増加する。𝐻e > 𝐻c2 では、内部の磁束密度は式(1.2)となり常伝導状態への相転移が起こる。内部の磁束密度がゼロの領域をマイスナー状態、磁束が量子化されて内部に侵入する領域を混合状態と呼ぶ。混合状態では電気抵抗がゼロに保たれるため、マイスナー状態と混合状態を合わせて超伝導状態と呼ぶ。

- 2 -

(a)

(b)

Fig. 1.1 外部磁界に対する超伝導体の内部の磁束密度

1.3Time-Dependent Ginzburg-Landau 方程式

ここでは、第 2 種超伝導体の磁気特性を記述する理論として、Ginzburg-Landau 理論[23]について述べる。

1.3.1 Ginzburg-Landau 方程式
超伝導状態とは電子が同位相である状態のことである。したがって量子力学的な性質
は巨視的なスケールまで保たれる。そこで量子力学の波動関数𝜓に対応する複素関数、
オーダーパラメータ𝛹を定義する。その大きさの 2 乗|𝛹 |2 は超伝導電子密度に比例する。
超伝導体の自由エネルギーを記述する。まず、超伝導電子密度に依存するエネルギー密
度𝐹𝛹 は、|𝛹 |2 で冪展開すると、
1
𝐹𝛹 = const. +𝛼 |𝛹 |2 + 𝛽 |𝛹 |2 + ⋯
2

(1.3)

となる。式(1.3)を凝縮エネルギー密度と呼ぶ。𝛼、𝛽は冪展開係数である。磁場を印加した場合、磁場のエネルギーも考慮する必要がある。磁場のエネルギー密度𝐹𝐵 は、
𝐹𝐵 =1(∇ × 𝑨)2
2𝜇0

(1.4)

と表される。ここで、𝑨はベクトルポテンシャルである。また、磁場印加時はオーダー
パラメータが空間的に変化することも考慮すべきである。これによって量子力学的な運動エネルギーも生じる。運動エネルギー密度𝐹𝐾 は、
𝐹𝐾 =1|(−iℏ∇ + 𝑒 ∗ 𝑨)𝛹 |2
2𝑚∗

(1.5)

と表される。ここでℏはディラック定数、𝑚∗ は超伝導電子の質量、𝑒 ∗ は超伝導電子の電荷である。これらのエネルギー密度、式(1.3)、(1.4)、(1.5)を足し合わせると磁場中の超伝導体の自由エネルギー密度が求まる。よって、超伝導体の自由エネルギー密度𝐹s (𝐵)

- 3 -

は、
1
1
1
𝐹s (𝐵) = 𝐹n (0) + 𝛼|𝛹|2 + 𝛽|𝛹|4 +(∇ × 𝑨)2 +|(−iℏ∇ + 𝑒 ∗ 𝑨)𝛹|2
2
2𝜇0
2𝑚∗

(1.6)

と記述できる。式(1.6)において、𝐹n (0)はゼロ磁界中の常伝導状態の自由エネルギー密度である。このとき、𝛹と𝑨は超伝導領域𝑉の全エネルギー𝐸𝑠 = ∫ 𝐹𝑠 d𝑉を最小にするように決定される。その方法として変分法を用いる。𝐸s を𝛹の共役複素数𝛹 ∗ と𝑨について変分し、その値が 0 となるようにすると、
δ𝐸s
∂𝐸s
𝜕𝐸s
=
−
[𝜵
∙
]=0
δ𝛹 ∗ ∂𝛹 ∗
𝜕𝜵𝛹 ∗

(1.7)

δ𝐸s 𝜕𝐸s
𝜕𝐸s=− [𝜵 ∙
]=0
δ𝑨
𝜕𝑨
𝜕𝜵𝑨

(1.8)

と書ける。(1.7)式、(1.8)式をそれぞれ解くと、
1(−iℏ∇ + 𝑒 ∗ 𝑨)2 𝛹 + 𝛼𝛹 + 𝛽 |𝛹 |2 𝛹 = 0
2𝑚∗
𝑱=

1
iℏ𝑒 ∗ ∗𝑒 ∗2∗) −|𝛹 |2 𝑨(𝛹∇×∇×𝑨=
∇𝛹
−
𝛹∇𝛹
𝜇0
2𝑚∗
𝑚∗

(1.9)
(1.10)となる。𝑱は電流密度である。式(1.9)、(1.10)を Ginzburg-Landau 方程式という。

1.3.2 Time-Dependent Ginzburg-Landau 方程式
GL 方程式に時定 数を導 入して時 間依 存性を持 たせた 方程式 は Time-Dependent
Ginzburg-Landau 方程式(以下 TDGL 方程式)[24]と呼ばれる。
𝛹と𝑨の時定数を𝜏Ψ 、𝜏A とおくと、(1.7)式と(1.8)式は、
δ𝐸s
∂𝐸s
𝜕𝐸s
𝜕𝜓
=
− [𝜵 ∙] = −𝜏Ψ
∗
∗
∗
δ𝛹
∂𝛹
𝜕𝜵𝛹
𝜕𝑡

(1.11)

δ𝐸s ∂𝐸s
𝜕𝐸s
𝜕𝑨=− [𝜵 ∙] = −𝜏A
δ𝑨
∂𝑨
𝜕𝜵𝑨
𝜕𝑡

(1.12)と書き換えることができる。(1.11)式、(1.12)式にゲージ変換∂𝜓 ∂𝜓
→+ i𝑒 ∗ 𝑉𝜓∂𝑡∂t(1.13)

∂𝑨 ∂𝑨
→
+ 𝜵𝑉
∂𝑡
∂𝑡

(1.14)を与える。𝑉はスカラーポテンシャルである。
𝜕𝛹1(−iħ𝜵 + 𝑒 ∗ 𝑨)2 𝛹 + 𝛼𝛹 + 𝛽 |𝛹 |2 𝛹 = 0𝜏Ψ (+ 𝑖𝑒 ∗ 𝑉𝛹) +
𝜕𝑡
2𝑚∗

- 4 -

(1.15)

𝜕𝑨1
iħ𝑒 ∗ ∗
𝑒 ∗2
∗) −
(𝛹
𝜏A ( + 𝜵𝑽) + 𝜵 × 𝜵 × 𝑨 +𝜵𝛹−𝛹𝜵𝛹|𝛹 |2 𝑨 = 0
𝜕𝑡
𝜇0
2𝑚∗
𝑚∗

(1.16)

超伝導体の数値解析は、式(1.15)と式(1.16)を解くことによって行うことができる。

1.4

磁束の量子化第 2 種超伝導体において、外部磁界が𝐻c1 < 𝐻e < 𝐻c2 のとき磁束は量子化された状態
で超伝導体の内部に入る。ここでは磁束の量子化について説明する。
オーダーパラメータを位相𝜑を用いて𝛹 = |𝛹 |i𝜑 と置くと、式(1.10)は
𝒋=−

2ℏ𝑒
4𝑒 2
2 𝛻𝜑 −|||𝛹 |2 𝑨𝛹𝑚∗𝑚∗

(1.17)と書ける。ここで、超伝導体内で孤立した単一の量子化磁束(中心𝑅 = 0)を考える。この磁束量子を取り囲み、中心から十分に離れた円C(𝑅 ≫ 𝜆)を考えると、この上では𝒋 = 0かつ|𝛹 |2 ≠ 0である。よって式(1.17)は、
𝑨=−

ℏ
𝛻𝜑
2𝑒

(1.18)と変形できる。C上で𝑨について線積分すると、これはCに鎖交する磁束𝛷に等しいため、∮ 𝑨 ∙ d𝒔 = ∮ 𝒃 ∙ d𝑺 = 𝛷C(1.19)

Sとなる。ここでSはCに囲まれた面である。式(1.19)に式(1.18)を代入すると、
𝛷=−ℏℏ∮ ∇𝜑 ∙ d𝒔 = − Δ𝜑2𝑒 C
2𝑒(1.20)となる。ただし、Δ𝜑はCを一周したときの位相の変化分である。𝛹は一価関数であるため、Δ𝜑は2𝜋の整数倍でなければならない。今は単一の量子化磁束の周りを考えている
ため、Δ𝜑 = −2𝜋、磁束は、
𝛷=ℎ= 𝜙0 = 2.0678 × 10−15 [Wb]2𝑒(1.21)となる。𝜙0 を磁束量子という[25]。

1.5

磁束ピンニング機構

超伝導現象には応用上で重大な問題がある。第 2 種超伝導体では、𝐻c1 以上の外部磁界が与えられると混合状態となり内部に量子化された磁束が侵入する。その状態の超伝- 5 -導体に電流を流すと、量子化磁束にローレンツ力𝑭L = 𝑱 × 𝑩が働き、磁束は速度𝒗で運動する。これにより誘導電界
𝑬= 𝑩×𝒗(1.22)が発生する。このため、電気抵抗が発生し、エネルギーが消費される。このままでは電気抵抗ゼロという超伝導の特長を活かすことができない。この問題は、磁束の運動を止めるピンを導入することで対策できる。磁束が止まれば𝒗 = 0より𝑬 = 0となるためである。ここでは、磁束のピンニング機構について説明する。

1.5.1 要素的ピン力
ピン止めの要因の一つに、常伝導析出物によるピンニングがある。常伝導析出物とは、ある元素を多量に含有させて超伝導体を作ることで超伝導体内にその元素が析出した
常伝導状態の部分のことである。
超伝導体を貫く磁束線の中心部(半径𝜉程度)は、ほとんど常伝導状態にあると言える。
つまり、周囲の超伝導部分に比べて凝縮エネルギーの分だけエネルギーは高い。超伝導体中に常伝導析出物があり、磁束線の常伝導核がそれと交差すると、交差しない場合よりも破壊される超伝導領域が減少し、凝縮エネルギーが、
1
𝑈p ≅ 𝜇0 𝐻c2 𝜋𝜉 2 𝐿
2(1.23)

だけ低くなる。ここで、𝐻c は熱力学的臨界磁界、𝜉はコヒーレンス長、1/2𝜇0 𝐻c2 は凝縮エ
ネルギー密度を表す。𝐿は、磁束線が常伝導析出物と交差する部分の高さである。その結果、磁束線と常伝導析出物との間に引力相互作用が働き、低いエネルギーを保つために磁束線は常伝導析出物の中で固定される。ここから要素的ピン力𝑓p (1 個のピンが及ぼす最大力)を計算すると、エネルギーの変化は2𝜉(常伝導核の直径)で起こるため、
𝑓p ≅

𝑈p 𝜋
≅ 𝜇 𝐻 2 𝜉𝐿
2𝜉 4 0 c(1.24)となる。

1.5.2 巨視的ピン力密度
高磁界下では、磁束線はたがいに接近しているため、相互にばねのような力が働き、1 本だけでなく複数の磁束線が一緒に移動する(磁束バンドル)。そこで、一辺の長さが𝐷の立方体で密度が𝑁p の常伝導析出物に複数の磁束線が交差する状態を考える。そのイメージ図を Fig. 1.2[25]に示す。

- 6 -Fig. 1.2 複数の磁束線のピン止め[25]
この図においてピンから力を受ける磁束線は赤で色を塗った 4 本のみである。磁束線
格子間隔を𝑎f とすると、ピン力を受ける磁束線の本数は𝐷/𝑎f であるので、この常伝導析出物によるピン力𝑓p は、
𝑓p =

𝜋
𝐷
𝜇0 𝐻c2 𝜉𝐷 ×
4
𝑎f

(1.25)と求められる。また、超伝導体全体のピン力密度(巨視的ピン力密度)𝐹p は、
𝜋𝜇0 𝐻c2 𝜉𝐷 2
𝐹p = 𝜂𝑁p 𝑓p = 𝜂𝑁p
4𝑎f(1.26)と求められる。ここで、𝜂はピンニング効率で、𝜂 < 1である。高磁界では、超伝導性が弱まるためオーダーパラメータの最大値が小さくなり、それにより凝縮エネルギーも小さくなる。その影響を式(1.26)に考慮すると、𝐹p は
𝐹p = 𝜂𝑁p

𝜋𝜇0 𝐻c2 𝜉𝐷 2
𝐵
(1 −
)
4𝑎f
𝜇0 𝐻c2(1.27)となる[25]。

1.6

Affine Integrator(AFI)

Affine Integrator[26, 27](以下、AFI)とは、拡散方程式やゲージ場存在下の Schrӧdinger
方程式、TDGL 方程式、回転系の Time-Dependent Gross-Pitaevskii(TDGP)方程式を数値積分するために考案された、陽的数値積分法である。
数値的安定性の点では、AFI はラプラシアン項に関して無条件安定であり、方程式全体に関して高い数値的安定性を持つ。また、Schrӧdinger 方程式のような線形保存系では
全エネルギーが厳密に保存され、TDGP 方程式のような非線形保存系にではノンドリフ
ト特性を持つ。
AFI は空間に関して離散化するための格子をチェッカボード状に分解することから導かれるアフィン変換対で構成される。数値的安定性とエネルギー保存の理論解析が容- 7 -易で、必要な記憶領域が最小限で済み、またワーキングメモリを追加する必要がないという利点を持つ。

1.6.1 1 次元系における空間・時間離散化
基礎方程式を 1 次元 TDGL 方程式
2
∂𝛹
𝜕
𝜏Ψ
= ( − i𝐴x ) 𝛹 − 𝛼𝛹 − 𝛽 |𝛹 |2 𝛹∂𝑡𝜕𝑥

(1.28)とする。ここではベクトルポテンシャルは定数与えられているとする。𝐴x はベクトルポ
テンシャルの x 成分である。
パラメータ𝜏𝛹 は時定数であり一般には複素数値をとる。𝜏𝛹 が実数値の場合は TDGL
方程式、𝜏𝛹 = −iのときは Time-Dependent Gross-Pitaevskii(TDGP)方程式あるいは非線形
Schrӧdinger 方程式と同じ形になる。
ベクトルポテンシャルはリンク変数[28] 𝑤𝑖𝑗 によって次のように実装される。
𝑤𝑖𝑗 = exp(i𝜃𝑖𝑗 ) ,𝑥𝑖 + 𝑥𝑗𝜃𝑖𝑗 = −ℎ𝐴x (),
2(1.29)

𝑖, 𝑗 = 1, 2, 3, 4

𝑥𝑖 は格子点 𝑖 の𝑥座標、ℎ は格子点間隔、𝜃𝑖𝑗 はリンク変数の偏角である。ただし、もし格子点 𝑖 と格子点 𝑗 が接続していない場合、 𝑤𝑖𝑗 = 0 であるとする。
空間に関する離散化のために導入した格子、リンク変数の状況を Fig. 1.3[29]に示す。

Fig. 1.3 格子・リンク変数[29](a)空間離散化格子(b)境界条件で更新される格子点とリンク変数

式(1.28)を空間に関して離散化する。 𝛹𝑖 (𝑡) = 𝛹(𝑥𝑖 , 𝑡),
けば、以下の連立微分方程式が得られる。

- 8 -𝛼𝑖 = 𝛼(𝑥𝑖 ),𝛽𝑖 = 𝛽(𝑥𝑖 )とお

𝜏Ψ

d𝛹1
1
= 2 (𝑤12 𝛹2 + 𝑤
̅41 𝛹4 − 2𝛹1 ) − 𝛼1 𝛹1 − 𝛽1 |𝛹1 |2 𝛹1d𝑡ℎ(1.30)𝜏Ψd𝛹3
1= 2 (𝑤34 𝛹4 + 𝑤
̅23 𝛹2 − 2𝛹3 ) − 𝛼3 𝛹3 − 𝛽2 |𝛹3 |2 𝛹3d𝑡ℎ(1.31)𝜏Ψd𝛹2
1= 2 (𝑤23 𝛹3 + 𝑤
̅12 𝛹1 − 2𝛹2 ) − 𝛼2 𝛹2 − 𝛽3 |𝛹2 |2 𝛹2d𝑡ℎ(1.32)𝜏Ψd𝛹41= 2 (𝑤41 𝛹1 + 𝑤
̅34 𝛹3 − 2𝛹4 ) − 𝛼4 𝛹4 − 𝛽4 |𝛹4 |2 𝛹4d𝑡ℎ(1.33)

ここで、𝑤
̅ は 𝑤 の複素共役を表す。ここで、
2

𝑈𝑖 = 𝛼𝑖 + 𝛽𝑖 |𝜓̃𝑖 | ,(1.34)

𝑖 = 1, 2, 3, 4

̃𝑖 は 𝛹𝑖 の推定値である。詳細は後述する。
とおく。𝛹
連立微分方程式は次のようにも書き表される。AFI 法の準備のため、奇数の格子点の組と偶数の格子点の組の順番に並べている。
−2 − 𝑈1 ℎ2
𝛹1
𝑑 𝛹
1
0
𝜏Ψ ( 3 ) = 2
𝑑𝑡 𝛹2
ℎ
𝑤
̅12
𝛹4
𝑤
̅14
(

0
−2 − 𝑈3 ℎ2
𝑤
̅32
𝑤
̅34

𝑤12
𝑤32
−2 − 𝑈2 ℎ2
0

𝑤14
𝑤34
0
−2 − 𝑈4 ℎ2 )(1.35)

式(1.35)は、次のようにも表される。
𝜏Ψ

ここで、𝒒 = (𝛹1𝑑 𝒒
1 𝐷(−𝜎𝑄 )(𝒑) = 2 (𝑑𝑡ℎ𝑊†𝛹3 )𝑇 、𝒑 = (𝛹2𝑊
𝐷(−𝜎𝑃 )𝒒) (𝒑)(1.36)

𝛹4 )𝑇 および

𝜎𝑖 = 2 + 𝑈𝑖 ℎ2 、

(1.37)

𝑖 = 1, 2, 3, 4

とおいた。また
𝐷(−𝜎Q ) = (

−𝜎1
0

0
−𝜎
) 、𝐷(−𝜎P ) = ( 2
−𝜎3
0

𝑤12
0
) 、𝑊 = (𝑤
−𝜎4
32

𝑤14
𝑤34 )

(1.38)

̅ 𝑇 である。
とおいた。なお、𝑊 † = 𝑊
式(1.35)あるいは式(1.36)は
𝛹1
𝛹1
𝑑 𝛹3
𝛹
𝜏Ψ ( ) = 𝐴̂ ( 3 )
𝛹2
𝑑𝑡 𝛹2
𝛹4
𝛹4

(1.39)あるいは
𝑑 𝒒𝒒( ) = 𝐴̂ (𝒑)𝑑𝑡 𝒑とあらわされる。ただし時間推進演算子 𝐴̂ は次のように定義される。
𝜏Ψ

- 9 -

(1.40)

𝐴̂ =

1
ℎ2

∑

(−𝜎𝑖 𝛹𝑖 +𝑖∈{1,2,3,4}∑

𝑤𝑖𝑗 𝜓𝑗 )𝑗∈{1,2,3,4}∂∂𝛹𝑖

時間推進演算子 𝐴̂ を次のように分解する。
𝐴̂ = 𝐴̂Q + 𝐴̂P(1.41)(1.42)

ただし、
𝐴̂Q =

1
𝜕
∑ (−𝜎𝑖 𝛹𝑖 + ∑ 𝑤𝑖𝑗 𝛹𝑗 )2ℎ𝜕𝛹𝑖(1.43)

1
𝜕
∑
(−𝜎
𝛹
+
∑
𝑤
𝛹
)
𝑖
𝑖
𝑖𝑗
𝑗
ℎ2
𝜕𝛹𝑖

(1.44)

𝑖∈𝑄

𝐴̂P =

𝑗∈𝑃

𝑖∈𝑃

𝑗∈𝑄

ここでは、𝑄 = {1, 3} および 𝑃 = {2, 4} である。
次に、微分方程式
𝜏Ψ

𝑑 𝒒
𝒒
𝒒
( ) = 𝐴̂ (𝒑) = (𝐴̂Q + 𝐴̂P ) (𝒑)𝑑𝑡 𝒑(1.45)

を時間に関して離散化する。時間刻み幅を𝜏とおくことで離散時間における時間発展方
程式
𝜏
𝜏
𝒒
𝒒
𝒒′
( ′ ) = exp ( 𝐴̂) (𝒑) = exp ( (𝐴̂Q + 𝐴̂P )) (𝒑)𝒑𝜏Ψ𝜏Ψ(1.46)

が得られる。なお、(𝒒, 𝒑) = (𝒒(𝑡), 𝒑(𝑡)) および(𝒒′ , 𝒑′ ) = (𝒒(𝑡 + 𝜏), 𝒑(𝑡 + 𝜏))である。
指数関数演算子は Lie-Trotter-Suzuki 分解により次のように近似できる。
𝜏𝜏𝐴̂𝑄 ) exp ( 𝐴̂𝑃 ) + 𝑂(𝜏 2 )
𝜏Ψ
𝜏Ψ
𝜏
𝜏
𝜏
𝜏
exp ( (𝐴̂𝑄 + 𝐴̂𝑃 )) = exp (𝐴̂𝑄 ) exp ( 𝐴̂𝑃 ) exp (𝐴̂ ) + 𝑂(𝜏 3 )
𝜏Ψ
2𝜏Ψ
𝜏Ψ
2𝜏Ψ 𝑄
⋯ + 𝑂(𝜏 4 )
(⋯
exp (

(1.47)

指数関数演算子 exp(𝜏Ψ −1 𝜏𝐴̂𝑄 ) および exp (𝜏Ψ −1 𝜏𝐴̂𝑃 ) の行列表現を得ることができれば数値積分スキーム(AFI)を構成できる。
まず、演算子 𝐴̂𝑄 を 𝛹𝑖 、𝑖 ∈ 𝑄 に作用させる。
𝐴̂𝑄 𝛹𝑖 =1(−𝜎𝑖 𝛹𝑖 + ∑ 𝑤𝑖𝑗 𝛹𝑗 ) ,ℎ2𝑗∈𝑃次に、𝐴̂𝑄 を 𝛹𝑖 、𝑖 ∈ 𝑄 に 2 回作用させる。

- 10 -

𝑖∈𝑄(1.48)

𝐴̂2𝑄 𝛹𝑖 = (−𝜎𝑖 2
1
) (𝛹𝑖 − ∑ 𝑤𝑖𝑗 𝛹𝑗 ) ,
2
ℎ
𝜎𝑖

𝑖∈𝑄

(1.49)𝑗∈𝑃

したがって、𝐴̂Q を 𝛹𝑖 、𝑖 ∈ 𝑄 に 𝑚回作用させると
𝐴̂𝑚𝑄 𝛹𝑖 = (−𝜎𝑖 𝑚1) (𝛹𝑖 − ∑ 𝑤𝑖𝑗 𝛹𝑗 ) ,2ℎ𝜎𝑖𝑖 ∈ 𝑄,

𝑚≥1(1.50)𝑗∈𝑃であることがわかる。ゆえに、
∞

𝑚
𝜏
1 𝜏
exp ( 𝐴̂𝑄 ) 𝛹𝑖 = 𝛹𝑖 + ∑( 𝐴̂𝑄 ) 𝛹𝑖
𝜏Ψ
𝑚! 𝜏Ψ
𝑚=1

= 𝑎𝑖 𝛹𝑖 + 𝑏𝑖 ∑ 𝑤𝑖𝑗 𝛹𝑗 ,𝑖∈𝑄

(1.51)

𝑗∈𝑃が得られる。ただし、
𝑎i = exp (−

𝜎𝑖 𝜏
),
ℎ2 𝜏Ψ

𝑏𝑖 =

1−𝑎𝑖
𝜎𝑖

(1.52)

とおいた。また、
𝐴̂𝑄 𝛹𝑖 = 0,(1.53)

𝑖∈𝑃

より、
𝐴̂𝑚
𝑖 ∈ 𝑃, 𝑚 ≥ 1
𝑄 𝛹 = 0,であることがわかる。演算子 𝐴̂P に関しても同様に、次の結果が得られる。

(1.54)

𝜏𝒒
𝐷(𝑎Q ) 𝐷(bQ )𝑊 𝒒
exp ( 𝐴̂𝑄 ) (𝒑) = () (𝒑 )
𝜏Ψ
0𝐼(1.55)

𝜏𝒒𝐼
exp ( 𝐴̂𝑃 ) (𝒑) = (
†
𝐷(𝑏
𝜏Ψ
P )𝑊

(1.56)𝒒0)( )𝐷(𝑎P) 𝒑

ここで、
𝑎
𝐷(𝑎𝑄 ) = ( 1
0

0
),
𝑎3

𝑏
𝐷(𝑏𝑄 ) = ( 1
0

0
)
𝑏3

𝑎
𝐷(𝑎𝑃 ) = ( 2
0

0
),
𝑎4

𝑏
𝐷(𝑏𝑃 ) = ( 2
0

0
)
𝑏4

(1.57)

とおいた。𝐼は単位行列である。
式(1.55)、式(1.56)はそれぞれ指数関数演算子 exp(𝜏𝛹 −1 𝜏𝐴̂𝑄 )およびexp(𝜏𝛹 −1 𝜏𝐴̂𝑃 )の行
列表現を表している。つまり、AFI が得られたことになる。

1.6.2 2 次元系での AFI
ここで、2 次元系を考える。基礎方程式は- 11 -𝜏Ψ∂𝜓= (∇ − i𝑨)2 𝛹 − 𝛼𝛹 − 𝛽|𝛹|2 𝛹∂𝑡(1.58)∂ ∂∇= ( , ) ,∂𝑥 ∂𝑦(1.59)と表される。ただし、𝑨 = (𝐴𝑥 , 𝐴y )である。ここでは格子点やリンク変数を指定するためにダブルインデックスを用いる。
このダブルインデックスは空間的な位置に対応する。Fig. 1.4[29]にダブルインデックスの定義を示す。リンク変数は𝑥成分と𝑦成分に分けて定義されることになる。

Fig. 1.4 2 次元系におけるダブルインデックスの定義(周期的境界条件)[29]
格子点は図中の黒丸と白丸のように、2 つのグループに分割される。まず空間に関して離散化すると次式のような微分方程式が得られる。
𝜏Ψ

d𝛹𝑖,𝑗
1
𝑦
𝑦
𝑥
𝑥
= 2 (𝑤𝑖,𝑗𝛹𝑖+1,𝑗 + 𝑤
̅𝑖−1,𝑗𝛹𝑖−1,𝑗 + 𝑤𝑖,𝑗 𝛹𝑖,𝑗+1 + 𝑤
̅𝑖,𝑗−1 𝛹𝑖,𝑗−1 − 4𝛹𝑖,𝑗 )d𝑡ℎ−𝛼𝑖,𝑗 𝛹𝑖,𝑗 − 𝛽𝑖,𝑗 |𝛹𝑖,𝑗 |2 𝛹𝑖,𝑗

(1.60)

ただし、周期的境界条件より、たとえば𝑖 = 0は𝑖 = 4、𝑖 = 5は𝑖 = 1と読み替えるものとする。ここでは、
̃𝑖,𝑗 |2
𝑈𝑖,𝑗 = 𝛼𝑖,𝑗 + 𝛽𝑖,𝑗 |𝛹(1.61)

とおく。ただし、𝜓̃𝑖,𝑗 =

1 𝑥
𝑦
𝑦
𝑥
(𝑤 𝛹
+𝑤
̅𝑖−1,𝑗𝛹𝑖−1,𝑗 + 𝑤𝑖,𝑗 𝛹𝑖,𝑗+1 + 𝑤
̅𝑖,𝑗−1 𝛹𝑖,𝑗−1 )4 𝑖,𝑗 𝑖+1,𝑗(1.62)である。また、
𝜎𝑖,𝑗 = 4 + 𝑈𝑖,𝑗 ℎ2- 12 -(1.63)

とおく。
1 次元系のときと同様、格子点を黒丸グループ𝑄と白丸グループ𝑃に分割することにより AFI を構成することができる。時間推進演算子𝐴̂ = 𝐴̂𝑄 + 𝐴̂𝑃 は 1 次元系のときと同様に定義することができる。黒丸グループに属する格子点上の𝛹の値で構成されるベク
トルを𝒒、白丸グループに属する格子点上の𝛹の値で構成されるベクトルを𝒑とおけば、
1 次元系のときと同様の次のような形の微分方程式
𝜏Ψ

d 𝒒
𝒒
𝒒
( ) = 𝐴 (𝒑) = (𝐴̂𝑄 + 𝐴̂𝑃 ) (𝒑)d𝑡 𝒑(1.64)

が得られる。この微分方程式を時間に関して時間刻み幅𝜏で離散化することで、離散時間における時間発展方程式
𝜏
𝜏
𝒒
𝒒
𝒒′
( ) = exp ( 𝐴̂) (𝒑) = exp ( (𝐴̂𝑄 + 𝐴̂𝑃 )) (𝒑)
𝒑′
𝜏Ψ
𝜏Ψ

(1.65)

が得られる。指数関数演算子の行列表現も 1 次元系の時と同様に、
exp (

𝜏
𝒒
𝐷(𝑎𝑄 ))( ) = (𝜏Ψ 𝒑0

𝐷(𝑏𝑄 )𝑊 𝒒
) (𝒑)
𝐼

𝜏
𝒒
𝐼
exp ( ) (𝒑) = (
†
𝐷(𝑏
𝜏Ψ
𝑃 )𝑊

𝒒
0
)( )
𝐷(𝑎𝑃 ) 𝒑(1.66)(1.67)

行列𝐷(𝑎𝑄 )および𝐷(𝑎𝑃 )は対角行列であり、それぞれ
𝑎𝑖,𝑗 = exp (−𝜎𝑖,𝑗 𝜏
𝜎𝑖,𝑗 𝜏
) , (𝑖, 𝑗) ∈ 𝑃 および 𝑎𝑖,𝑗 = exp (− 2 ) , (𝑖, 𝑗) ∈ 𝑃
2
ℎ 𝜏Ψ
ℎ 𝜏Ψ

(1.68)を対角成分とするものである。なお、ダブルインデックスにより指定される格子点が黒丸グループに属することを(𝑖, 𝑗) ∈ 𝑃、白丸グループに属することを(𝑖, 𝑗) ∈ 𝑄と表記している。
また、行列𝐷(𝑏𝑄 )および𝐷(𝑏𝑃 )は対角行列であり、それぞれ
𝑏𝑖,𝑗 =

1 − 𝑎𝑖,𝑗
1 − 𝑎𝑖,𝑗, (𝑖, 𝑗) ∈ 𝑄 および 𝑏𝑖,𝑗 =, (𝑖, 𝑗) ∈ 𝑃
𝜎𝑖,𝑗
𝜎𝑖,𝑗

(1.69)を対角成分とするものである。
行列𝑊は黒丸グループに属する格子点と白丸グループに属する格子点を結ぶ辺の上で定義されたリンク変数を成分とする行列である。行列𝑊の行は白丸グループの格子点に対応し列は黒丸グループの格子点に対応する[5]。

- 13 -

1.7

リザバーコンピューティング

リザバーコンピューティングは、特に系列データや時系列データの機械学習に適した計算フレームワークである。その特徴として、学習に必要な計算量が少ないため、よく知られているディープラーニングよりも高速な計算が可能であること、純粋に通常の計算機で扱えることなどの利点がある。また、リザバーコンピューティングは、デジタル計算機だけでなく、様々なスケールの物理システムを用いたハードウェア上でも実装可能であるという側面からも、近年注目されている[30]。

1.7.1 リザバーコンピューティングの概念
リザバーコンピューティングは、いくつかの特定のリカレントニューラルネットワー
クから派生した一般的な概念である。リザバーコンピューティングモデルは一般的に、リザバー（Reservoir）とリードアウト（Readout）から構成される。リザバーは固定されており、学習プロセス中に変化することはない。一方、リードアウトは線形学習機械のような単純な学習識別器を使用することで、高速な学習を実現している。。
リザバーコンピューティングの概念図を Fig. 1.5 に示す。英語の "reservoir "は液体を入れる容器、または日常使用する水を貯める貯水池を意味する。例えば、水面（リザバー）に石を落とすと波紋が広がる。この波紋は、投入された石の大きさや形によって変化する。つまり、水面の波紋は入力された情報を反映していると考えられる。
ここで、最初の石を投入して波紋が残っている間に、別の石を次々に投入することを考える。ここで、順番に石を投げる順序が時系列の入力に対応する。複数の石が入力されると、それぞれの石から発生する波紋は互いに影響し合い、複雑なものとなる。波紋は石の大きさや形だけでなく、石を入れる順番にも影響される。これは、波紋がその時点までに入力された石に関する情報、つまり過去の入力に関する情報を「記憶」していることを意味する。この性質が、時系列データを扱うための鍵となる。リザバーコンピューティングのアイデアは、入力された石の情報を波紋から読み出し、時系列予測を行うことである。

Fig. 1.6 リザバーコンピューティングの概念図- 14 -

1.7.2 エコーステートネットワーク
エコーステートネットワーク（ESN）はリザーバーコンピューティングの代表的なモデルの一つである。これは、人工ニューラルネットワークをベースとし、固定結合重みを持つリカレント・ニューラル・ネットワーク（reservoir）を使用して、時系列入力の過去情報がエコーされた状態（Echo state）を作り出し、そこから入力の特徴を読み出す（read out）というものである。
エコーステートネットワークの最も基本的なモデルを以下に説明する。モデルの構造を Fig. 1.6 に示す。ここで、𝑊 in 、𝑊、𝑊 out はそれぞれ入力層からリザバー層への結合
重み、リザバー内のノード間の結合重み、リザバー層から出力層への結合重みである。
一般的なリカレントニューラルネットワークでは𝑊in 、𝑊を学習させる必要であるが、エコーステートネットワークでは𝑊in 、𝑊は固定である。したがって𝑊out だけが学習ア
ルゴリズムによって調整される。リザバーは入力データの変換器であり、リードアウト
はリザバーの状態から入力の特徴を適切に読み出すための学習器である。

Fig. 1.7 エコーステートネットワークの基本モデル
Fig. 1.6 において、入力層のノード数を𝑁𝑢 、リザバーのノード数を𝑁𝑥 、出力層のノード数を𝑁𝑦 と表す。時系列データを扱うため、各ノードの状態は離散時間𝑛(𝑛 = 0, 1, 2 … )とともに発展すると仮定する。時刻𝑛における入力ベクトル𝒖(𝑛)、リザバーのノード状
態のベクトル𝒙(𝑛)、出力ベクトル𝒚(𝑛)をそれぞれ、⊤

𝒖(𝑛) = (𝑢1 (𝑛), … , 𝑢𝑁𝑢 (𝑛)) ∈ ℝ𝑁𝑢⊤𝒙(𝑛) = (𝑥1 (𝑛), … , 𝑥𝑁𝑢 (𝑛)) ∈ ℝ𝑁𝑢⊤𝒚(𝑛) = (𝑦1 (𝑛), … , 𝑦𝑁𝑢 (𝑛)) ∈ ℝ𝑁𝑢(1.70)(1.71)(1.72)と表す。また、それぞれの結合重み行列を、
in
𝑊in = (𝑤𝑖𝑗) ∈ ℝ𝑁𝑥 ×𝑁𝑢(1.73)) ∈ ℝ𝑁𝑥 ×𝑁𝑥(1.74)

𝑊 = (𝑤𝑖𝑗- 15 -

𝑊out = (𝑤𝑖𝑗out) ∈ ℝ𝑁𝑥 ×𝑁𝑦(1.75)と表す。
各時刻𝑛において、入力ベクトルが入力層に与えられると、その情報は変換されながら、出力層に向かって伝達される(順伝播)。時刻𝑛のリザバー状態ベクトル𝒙(𝑛)と時刻(𝑛 + 1)の入力ベクトル𝒖(𝑛 + 1)が与えられると、リザバー状態ベクトルの時間発展は、
𝒙(𝑛 + 1) = 𝒇 (𝑊 in 𝒖(𝑛 + 1) + 𝑊𝒙(𝑛)) (𝑛 = 0, 1, 2, … )(1.76)と記述される。これをリザバーのノード𝑖(𝑖 = 1, … , 𝑁𝑥 )について書き下すと、
𝑁𝑢𝑁𝑥𝑥𝑖 (𝑛 + 1) = 𝑓 (∑ 𝑤𝑖𝑗𝑖𝑛 𝑢𝑗 (𝑛 + 1) + ∑ 𝑤𝑖𝑗 𝑥𝑗 (𝑛))𝑗=1
𝑗=1

(1.77)となる。ここで、𝒇(∙)は活性化関数𝑓を括弧内の要素ごとに施す操作を表す。出力層の出力ベクトルは、リザバーのすべてのノード状態の線形結合として、
𝒚(𝑛 + 1) = 𝑊 out 𝒙(𝑛 + 1) (𝑛 = 0, 1, 2, … )(1.78)と与えられる。これを出力ノード𝑘(𝑘 = 1, … , 𝑁𝑦 )について書き下すと、
𝑁𝑥out (𝑛𝑦𝑘 (𝑛 + 1) = ∑ 𝑤𝑘𝑖
𝑥𝑖 + 1)(1.79)

𝑗=1となる。リザバー状態ベクトルの初期条件𝒙(0)と入力ベクトル𝒖(1)が与えられると、式(1.76)および式(1.78)から、リザバー状態ベクトル𝒙(1)と出力ベクトル𝒚(1)が計算される。同様の状態更新により、リザバー状態ベクトルの時間発展𝒙(𝑛)とモデル出力の時系
列𝒚(𝑛)が生成される。

1.7.3 バッチ学習(線形回帰とリッジ回帰)
バッチ学習では、データを一定期間貯めておき、それらを用いて学習を一括で行う。
その時間範囲内では𝑊 out は時間𝑛によらず一定とする。ここでは線形回帰とリッジ回帰について説明する。
エコーステートネットワークの基本モデルの出力は、式より、𝒚(𝑛) = 𝑊 out 𝒙(𝑛)と表される。理想的には、すべての時刻𝑛でこれが目標出力𝒅(𝑛)と一致すればよく、𝑊 out 𝒙(𝑛) = 𝒅(𝑛) (𝑛 = 1, … , 𝑇)

(1.80)

を満たす𝑊 out ∈ ℝ𝑁𝑦 ×𝑁𝑥 を求めればよい。これは、連立線型方程式を解くことに他ならない。
ここで、時刻𝑛 = 1, … , 𝑇についてリザバー状態ベクトルを横方向に連結した行列を、
𝑋 = [𝒙(1), … , 𝒙(𝑇)] ∈ ℝ𝑁𝑥 ×𝑇(1.81)と定義し、目標時系列出力ベクトルを連結した行列を、
𝐷 = [𝒅(1), … , 𝒅(𝑇)] ∈ ℝ𝑁𝑥 ×𝑇- 16 -(1.82)

と定義すると、式(1.80)は、
𝑊 out 𝑋 = 𝐷(1.83)
̂ out =と書ける。もし𝑁𝑥 = 𝑇で𝑋が正則行列ならば、その逆行列𝑋 −1 を用いて解を 𝑊
𝐷𝑋 −1 と求めることができる。
しかし、現実の多くの問題では𝑇 ≫ 𝑁𝑥 であり、式(1.83)は未知変数の数より制約式の数が多い優決定系となるので、近似解を求める必要がある。そこで、式の両辺を最小化するため、線形回帰を行う。具体的には、最小二乗法を利用して、二乗誤差の総和
1
2
𝐸LR = ||𝐷 − 𝑊 out 𝑋||𝐹2
𝑇1
2
= ∑||𝒅(𝑛) − 𝒚(𝑛)||2
2
𝑛=1
𝑇

1
2
= ∑||𝒅(𝑛) − 𝑊 out 𝒙(𝑛)||22(1.84)

𝑛=1を最小化する。ここで、|| ∙ ||𝐹 は行列の Frobenius ノルムを表す。この最小化問題の解は、
正規方程式
𝑊 out 𝑋𝑋⊤ = 𝐷𝑋⊤(1.85)

̂ out = 𝐷𝑋⊤ (𝑋𝑋⊤ )−1𝑊(1.86)

を解くことにより、と求められる。したがって、𝑋の Moore-Penrose 疑似逆行列𝑋 † = 𝑋⊤ (𝑋𝑋⊤ )−1 を用いて̂ out = 𝐷𝑋 †(1.87)
𝑊と書ける。𝑋のサイズが大きい場合は、𝑋 † を直接求めると計算に必要なメモリ量が増加
して効率的ではないため、𝑋𝑋⊤ と𝐷𝑋 ⊤ を𝑛の増加とともに逐次的に求めて式(1.86)を適用するとよい。𝑋𝑋⊤ が正則でない場合には̂ out = 𝐷𝑋⊤ (𝑋𝑋 ⊤ )†𝑊

(1.88)

としてもよい。
リードアウトの学習パラメータ数はリザバーのノード数𝑁𝑥 に比例するため、𝑁𝑥 を増やせば出力誤差はより小さくなると予想される。しかし、モデルの自由度が高すぎると、モデルが訓練データに過剰適応してしまう(過学習になる)可能性がある。機械学習では、過学習を抑制するために正則化がよく用いられる。正則化は出力重み行列の要素の絶対値を小さくする傾向があるため、ノイズに強いモデルを得ることが期待できる。このよ
うに、必要最小限の学習パラメータでモデルを表現しようとする手法をスパース最適化と呼ぶ。
式(1.84)の右辺に学習パラメータの二乗和を正則化項として加える場合をL2 最適化と呼び、この場合の最適な𝑊 out を求める回帰をリッジ回帰という。すなわち、最小化す- 17 -

べきコスト関数を、
𝐸RR =

1
𝛽
2
2
||𝐷 − 𝑊 out 𝑋||𝐹 + ||𝑊 out ||𝐹
2
2
𝑇

𝑁𝑦 𝑁𝑥

𝑛=1

𝑖=1 𝑗=11𝛽22= ∑||𝒅(𝑛) − 𝑊 out 𝒙(𝑛)||2 + ∑ ∑|𝑤𝑖𝑗 |2
2(1.89)とする。ここで、𝛽 > 0は正則化項の大きさを調整する正則化パラメータである。このとき解は、
̂ out = 𝐷𝑋⊤ (𝑋𝑋⊤ + 𝛽𝐼)−1𝑊(1.90)と求められる。ここで、𝐼 ∈ ℝ𝑁𝑥 ×𝑁𝑥 は単位行列を表す。𝛽 = 0のときは式(1.86)に帰着される。

1.7.4 評価指標
̂ out をリードアウト
訓練データと学習アルゴリズムから求められた出力結合重み行列𝑊に用いるモデルを学習済みのモデルと呼ぶ。学習済みのモデルの性能は、出力 𝑦(𝑛)と目
標出力データ𝑑(𝑛)との間の誤差または一致度で評価される。ここでは Normalized Mean
Square Error(NMSE)と R Squared(R2)の 2 つを紹介する。NMSE と R2 はそれぞれ、
NMSE =∑𝑇𝑛=1( 𝑦(𝑛) − 𝑑(𝑛))2∑𝑛𝑘=1(𝑑(𝑛))2(1.91)∑𝑇𝑛=1(𝑦(𝑛) − 𝑑(𝑛))2∑𝑇𝑛=1(𝑑(𝑛) − 𝑌̅)2

(1.92)

R2 = 1 −と定義される。ここで 𝑌̅は出力の平均である。NMSE は誤差、R2 は一致度を表す指標である。

1.7.5 NARMA モデル
時系列予測とは、時系列データの将来値を予測する問題である。ここでは、非線形自己回帰移動平均(Nonlinear Auto Regressive Moving Average, NARMA)モデル[31]について説明する。このモデルは、自己回帰部分と移動平均部分からなる自己回帰移動平均（ARMA）モデルの一種であり、現在値と過去値の間に非線形依存性がある。
時系列入力𝑢(𝑛)が与えられるとき、時数𝑚の NARMA モデルの出力𝑑(𝑛 + 1)は、
𝑚−1

𝑑(𝑛 + 1) = 𝑎1 𝑑(𝑛) + 𝑎2 𝑑(𝑛) ∑ 𝑑(𝑛 − 𝑖) + 𝑎3 𝑢(𝑛 − 𝑚 + 1)𝑢(𝑛) + 𝑎4

(1.93)

𝑖=0と記述される。ここで、𝑚は 1 以上の整数、𝑎1 , 𝑎2 , 𝑎3 , 𝑎4 は定数パラメータである。右辺の第 2 項は過去𝑚ステップの状態に依存する非線形項である。右辺の第三項は入力に依存する非線形項である。NARMA モデルの時系列予測を行うためには、機械学習モデル- 18 -

が非線形性を持ち、少なくとも𝑚ステップ過去の入力を「記憶」を記憶している必要がある。𝑚 = 2のときのモデルを NARMA2、𝑚 = 10のときのモデルを NARMA10 という。

1.7.6 非線形-メモリタスク
リザバーコンピューティングのタスクの 1 つに非線形-メモリタスク[32]と呼ばれるものがある。このタスクの教師データは以下の目的関数、
𝑦(𝑛) = sin(𝜈 × 𝑢(𝑛 − 𝜏))(1.94)

によって生成される。ここで𝜈は非線形性をあらわるパラメータであり、𝜏は記憶性を表すパラメータである。𝜈と𝜏を変化させながら誤差または一致度を測定することで、リザバーの持つ非線形性と記憶性を同時に測ることができる。

1.8

研究目的

ここで本研究の目的を提示する。本研究の目的は主に 2 つある。第 2 種超伝導体の𝐸-𝐽特性の非線形性がリザバーコンピューティングに適しているということは卒業論文で示された。しかし、超伝導リザバーは実際に時系列予測を行うには十分な精度を得られていないということが課題として残っている。非線形性や記憶性に関しても改善の余地が見られた。したがって 1 つ目の目的は、卒業研究に用いたリザバーコンピュ
ーティングのタスクと同様のタスクを超伝導体側のパラメータや条件を変えながら実行し、精度を測定していくことで、超伝導リザバーの予測精度を向上させることである。また、1.1 節でも述べたように超伝導リザバーに限らず、リザバーコンピューティ
ングに全体の課題としてリザバー内部がブラックボックスであるため、内部の動きと時系列予測の結果との関係が不明瞭で、より高精度のリザバーを設計する確かな方法が確立されていないというものがある。本研究に用いる 2 次元 TDGL 方程式の数値シ
ミュレーションは超伝導体内部の磁束線の運動、それに伴う電磁現象、ピンの並びなどを実際に見ながら行うことができるという利点がある。つまり、リザバーの内部状
態と時系列予測の結果を同時に観察することができる。したがって 2 つ目の目的は、各パラメータを変更したことによる超伝導体内部に生じる変化を観察し、それに伴い予測精度がどのように変化するのかを明らかにすることで、超伝導リザバーにおけるリザバー内部と時系列予測の結果を関連付けることである。予測精度向上とリザバーと結果の関係の明瞭化が達成されれば、超伝導リザバーの他のリザバーに対する優位性を示すことにつながる。また、リザバーコンピューティング全体の持つ課題解決の糸口にもづながることを期待する。

- 19 -

第2章
2.1

実装および計算方法

実装する TDGL 方程式と固定パラメータの規格化

2 次元 TDGL 方程式を AFI 法によって実装した。本研究ではスカラーポテンシャルが常にいたるところでゼロであるゲージを採用した。ゲージ場(ベクトルポテンシャル)が
時間的、空間的ともに変化するとき、規格化 TDGL 方程式は以下の 2 つの方程式、
𝜕𝛹= (∇ − i𝑨)2 𝛹 − 𝛼𝛹 − 𝛽|𝛹 |2 𝛹
𝜕𝑡

(2.1)

𝜕𝑨
̅ (∇ − i𝑨)𝛹] − ∇ × ∇ × 𝑨= 𝐼𝑚[𝛹
𝜕𝑡

(2.2)

𝜏Ψ
𝜏A

で与えられる。式(2.1)がオーダーパラメータ、式(2.2)がベクトルポテンシャルに関する方程式である。ここで式(2.2)の𝜏A はベクトルポテンシャルに関する時定数である。冪展開係数𝛼と𝛽はそれぞれ、𝛼 = −20、𝛽 = 20と設定した。これによりコヒーレンス長は、1
𝜉=√≅ 0.22|𝛼 |(2.3)に規格化される。長さのパラメータをすべて𝜉で規格化した。したがって磁場侵入長と空間刻み幅は𝜉を使ってそれぞれ、
1𝛽𝜆=√= −√ = 1.0 ≅ 4.5𝜉
2
|𝛹∞ |
𝛼

(2.4)

ℎ = 0.1𝜉(2.5)と規格化された。また、規格化 TDGL 方程式において、上部臨界磁束密度𝐵c2 は、
𝐵c2 = √2|𝛼 |(2.6)と規格化した。また、シミュレーション領域はlength × width = 10ξ × 10ξとした。印加
磁場は計算領域に対して垂直に𝐵 = 0.035𝐵c2 で与えた。

- 20 -

2.2

空間離散化

ゲージ場の時間発展に関しては、ベクトルポテンシャル𝑨の形のまま実装するのでは
𝑖,𝑗

𝑖,𝑗

なく、容易にするためにリンク変数𝑤𝑥 , 𝑤𝑦 に関する偏微分方程式を導き、それを時間に関して離散化する方法を採用した。

Fig. 2.1 演算∇ ×のリンク変数による差分化およびループ変数の定義[29]リンク変数の積をループ変数𝐶と呼ぶこととする。このとき、式(2.2)における演算∇ × ∇ ×の差分化はループ変数𝐶の積で実現できる。
𝑖,𝑗

𝑖,𝑗

ここで、リンク変数𝑤x , 𝑤y の更新は次のように行った。
方程式(2.2)を空間に関して離散化することを考える。リンク変数の定義式
𝑖,𝑗

𝑖,𝑗

𝑤𝑥𝑖 ≡ exp(−iℎ𝐴𝑥 ) ,

𝑖,𝑗

(2.7)

𝑤𝑦 ≡ exp (−iℎ𝐴𝑦 )および Fig. 2.1(a)より、格子点(𝑖, 𝑗)を左下の頂点とする微小正方形によって演算子∇ ×を𝑖,𝑗

差分化することを考える。ループ変数𝐶z を次式のように、微小正方形の辺上で定義されたリンク変数の積で定義する。
𝑖,𝑗

𝑖+1,𝑗

𝐶𝑧 ≡ 𝑤𝑥 𝑤𝑦

𝑖.𝑗+1

𝑖,𝑗

𝑤
̅𝑥

𝑤
̅𝑦

𝑖+1,𝑗

− 𝐴𝑦 − 𝐴𝑥= exp[−iℎ(𝐴𝑦

𝑖,𝑗

𝑖,𝑗

≅ exp [−iℎ2 (𝛻 × 𝑨)𝑧 ]
𝑖,𝑗

1
2

1
2

𝑖,𝑗+1

𝑖,𝑗

+ 𝐴𝑥 )](2.8)

ここで、(∇ × 𝑨)𝑥 は、位置(𝑖 + , 𝑗 + )における𝑨のローテーションの𝑧成分である。

- 21 -

𝑖,𝑗

さらに、このループ変数𝐶z を用いて次式が得られる。
𝑖,𝑗 𝑖,𝑗−1
𝑖,𝑗𝐶z 𝐶z̅
≅ exp[−iℎ3 (∇ × ∇ × 𝑨)𝑥 ](2.9)

𝑖,𝑗 𝑖−1,𝑗𝑖,𝑗(2.10)𝐶z̅ 𝐶z≅ exp[−iℎ3 (∇ × ∇ × 𝑨)𝑦 ]𝑖,𝑗𝑖,𝑗
ここで式(2.9)と(2.10)は、それぞれリンク変数𝑤x および𝑤y と同じ位置で定義されてい

る。
𝑖,𝑗

𝑖,𝑗

𝑖,𝑗

𝑖,𝑗−1

ループ変数𝐶𝑧 は(∇ × 𝑨)𝑧 に対応する(Fig. 2.1[5])。Fig. 2.1(b)、Fig. 2.1(c)は𝐶𝑧 𝐶𝑧̅

と𝑖,𝑗 𝑖−1,𝑗𝐶𝑧̅ 𝐶𝑧の定義を説明している。

以上の計算より、TDGL 方程式を空間に関して離散化することにより、次の常微分方程式が得られた。
γ

∂𝛹𝑖,𝑗
1
𝑖,𝑗
𝑖−1,𝑗
𝑖,𝑗
𝑖,𝑗−1
= 2 (𝑤𝑥 𝛹𝑖+1,𝑗 + 𝑤
̅𝑥 𝛹𝑖−1,𝑗 + 𝑤𝑦 𝛹𝑖,𝑗+1 + 𝑤
̅𝑦 𝛹𝑖,𝑗−1 − 4𝛹𝑖,𝑗 )∂𝑡ℎ̃𝑖,𝑗 |2 𝛹𝑖,𝑗
−𝛼𝑖,𝑗 𝛹𝑖,𝑗 − 𝛽𝑖,𝑗 |𝛹(2.11)

𝑖,𝑗∂𝑤𝑥
1
̅𝑖,𝑗 𝑤𝑥𝑖,𝑗 𝛹𝑖+1,𝑗 ]𝑤𝑥𝑖,𝑗 − (𝐶𝑧̅ 𝑖,𝑗−1 𝐶𝑧𝑖,𝑗 − 1)𝑤𝑥𝑖,𝑗𝜏A= −iIm[𝛹
∂𝑡
ℎ2
𝑖,𝑗

𝜏A

∂𝑤𝑦
∂𝑡

̅𝑖,𝑗 𝑤𝑦𝑖,𝑗 𝛹𝑖+1,𝑗 ]𝑤𝑦𝑖,𝑗 −= −iIm[𝛹1𝑖,𝑗 𝑖−1,𝑗𝑖,𝑗(𝐶̅ 𝐶− 1)𝑤𝑦ℎ2 𝑧 𝑧(2.12)(2.13)

ここで、(∇ × ∇ × 𝑨)𝑖,𝑗𝑥 ≅𝑖,𝑗−1 𝑖,𝑗
𝐶𝑧 − 1
−iℎ3

𝐶𝑧̅

𝑖,𝑗

(∇ × ∇ × 𝑨)𝑖,𝑗
𝑦 ≅

𝑖−1,𝑗

𝐶𝑧̅ 𝐶𝑧−13−iℎであることを用いた。

- 22 -(2.14)(2.15)

2.3

境界条件の処理

ここでは境界条件の実装方法について説明する。
これまでのように 2 次元では長方形シミュレーション領域を考える。オーダーパラメー
タの値は格子点上で定義されている：
𝛹𝑖,𝑗 ,

𝑖 = 0, ⋯ , 𝑁𝑥 + 1,

𝑗 = 0, ⋯ , 𝑁𝑦 + 1(2.16)

ここでは𝑧方向に印加されている一様磁界𝐵a を考える。このとき、∇ × 𝑨 = 𝑩a =(0, 0, 𝐵𝑧 )が東西南北の境界で満たされていなくてはならない。さらに、𝑦方向に電流𝑱a が印加されている場合、すなわち𝑱a = (0, 𝐽a , 0)である場合、先の境界条件に電流による効
果が重畳されることになる。Fig. 2.2[29]に境界付近の格子点およびリンク変数の配置を示す。印加磁界および印加電流はリンク変数の境界条件に反映される。

Fig. 2.2 境界付近の格子点およびリンク変数の配置[29](a)西側境界(𝑥 = 0)(b)東側境界(𝑥 = 𝐿𝑥 ) (c)南側境界(𝑦 = 0) (d)北側境界(y = 𝐿𝑦 )- 23 -

2.3.1 リンク変数の境界条件
西側および東側の境界(𝑥 = 0および𝑥 = 𝐿𝑥 )において、境界面に垂直な成分𝑤xはノイマ
ン境界条件
0,𝑗

1,𝑗

𝑤𝑥 𝐵 = 𝑤𝑥 𝐵 ,(2.17)𝑁 ,𝑗
𝑁 −1,𝑗𝑤𝑥 𝑥 𝐵 = 𝑤𝑥 𝑥 𝐵(2.18)を採用した。ここで、𝑗B = 1, ⋯ 𝑁𝑦 である。右辺が左辺に代入されるような実装を行うものとし、以降の境界条件の式についても、同様に実装を行った。
西側の境界(𝑥 = 0)において境界面に平行な成分𝑤𝑦 が満たすべき条件は、アンペールの法則より、0,𝑗1,𝑗0,𝑗 +1

𝑤𝑥 𝐵 𝑤𝑦 𝐵 𝑤
̅𝑥 𝐵

0,𝑗

𝑤
̅𝑦 𝐵 ≅ exp (−i ∮𝑨 ∙ d𝒔)(0,𝑗𝐵)= exp (−i ∫ ∫∇ × 𝑨 ∙ d𝑺)(0,𝑗𝐵 )

≅ exp [−iℎ2 (𝐵𝑎 + 𝜇0

𝐿𝑥
𝐽 )]
2 𝑎

(2.19)

ここで𝜇0 は真空透磁率、𝐿𝑥 はシミュレーション領域の𝑥方向のサイズに対応する。
式(2.19)により、西側境界における境界面に平行な成分𝑤𝑦 として、
0,𝑗

0,𝑗

1,𝑗

0,𝑗 +1

𝑤𝑦 𝐵 = 𝑤𝑥 𝐵 𝑤𝑦 𝐵 𝑤
̅𝑥 𝐵

exp [iℎ2 (𝐵𝑎 + 𝜇0

𝐿𝑥
𝐽 )]
2 𝑎

(2.20)が得られた。
東側の境界(𝑥 = 𝐿𝑥 )において境界面に平行な成分𝑤𝑦 が満たすべき条件は、同様にアンペールの法則から導かれる次の関係式によって求めた。
𝑁 ,𝑗

𝑁 +1,𝑗𝐵

𝑤𝑥 𝑥 𝐵 𝑤𝑦 𝑥

𝑁 ,𝑗 +1𝑤
̅𝑥 𝑥 𝐵

𝑁 ,𝑗𝑤
̅𝑦 𝑥 𝐵 ≅ exp (−i ∮𝑨 ∙ d𝒔)(𝑁𝑥 ,𝑗𝐵)= exp (−i ∫ ∫∇ × 𝑨 ∙ d𝑺)(𝑁𝑥 ,𝑗𝐵)≅ exp [−iℎ2 (𝐵𝑎 − 𝜇0

𝐿𝑥
𝐽 )]
2 𝑎

(2.21)

式(2.21)により、東側境界における境界面に平行な成分𝑤𝑦 として、
𝑁 +1,𝑗𝐵

𝑤𝑦 𝑥

𝑁 ,𝑗

𝑁 ,𝑗 +1= 𝑤𝑥 𝑥 𝐵 𝑤𝑥 𝑥 𝐵

𝑁 ,𝑗𝑤𝑦 𝑥 𝐵 exp [−iℎ2 (𝐵𝑎 − 𝜇0

𝐿𝑥
𝐽 )]
2 𝑎

(2.22)が得られた。
ここで、式(2.20)および式(2.22)において、𝑗𝐵 = 1, ⋯ 𝑁𝑦 − 1である。
印加磁界(𝑧方向)と印加電流(𝑦方向)の効果は東西境界面における面に平行な成分に反映される。印加磁界と印加電流のリンク変数の境界条件への反映について、Fig. 2.3[29]- 24 -

に示す。

Fig. 2.3 東西境界における印加磁界、印加電流のリンク変数への反映[29]
南側、北側の境界面における境界面に垂直な成分に関して、次のようなノイマン境界条
件
𝑖 ,0

𝑖 ,1

𝑤𝑦𝐵 = 𝑤𝑦𝐵
𝑖𝐵,𝑁𝑦

𝑤𝑦

𝑖𝐵,𝑁𝑦 −1= 𝑤𝑦

(2.23)(2.24)を採用した。ここで、𝑖𝐵 = 1, ⋯ , 𝑁𝑥 である。
南および北側の境界面における接線成分についてもノイマン境界条件が満たされる
とした。
𝑖 ,0

𝑖 ,1

𝑤𝑥𝐵 = 𝑤𝑥𝐵

(2.25)

𝑖𝐵,𝑁𝑦 +1
𝑖𝐵,𝑁𝑦
𝑤𝑥
= 𝑤𝑥

(2.26)

ここで𝑖𝐵 = 1, ⋯ , 𝑁𝑥 − 1である。

- 25 -

2.3.2 オーダーパラメータの境界条件
西側と東側の境界においてはノイマン境界条件に対応する、次式を実装した。
0,𝑗

𝛹0,𝑗𝐵 = 𝑤𝑥 𝐵 𝛹1,𝑗𝐵(2.27)𝑁 ,𝑗
𝛹𝑁𝑥 +1,𝑗𝐵 = 𝑤
̅𝑥 𝑥 𝐵 𝛹𝑁𝑥,𝑗𝐵(2.28)

ここで、𝑗𝐵 = 1, ⋯ , 𝑁𝑦 である。
南側と北側の境界においても、ノイマン境界条件に対応する、次式を実装した。
𝑖 ,0

𝛹𝑖𝐵 ,0 = 𝑤𝑦𝐵 𝛹𝑖𝐵 ,1𝑖𝐵,𝑁𝑦𝛹𝑖𝐵 ,𝑁𝑦+1 = 𝑤
̅𝑦𝛹𝑖𝐵,𝑁𝑦(2.29)(2.30)

ここで、𝑖𝐵 = 1, ⋯ , 𝑁𝑥 である。
オーダーパラメータについて、値の更新方法を示す図を Fig. 2.4[29]に示す。

Fig. 2.4 境界条件と基礎方程式によるオーダーパラメータの更新[29]曲線で囲まれている格子点およびリンク変数は境界条件によって、囲まれていないもの
は基礎方程式によって値が更新される。

- 26 -

2.4

初期条件と時間離散化

オーダーパラメータに関する初期条件として、
𝛹(𝑥, 𝑦, 0) = cos (𝜋𝑚𝜋𝑛𝑥) cos ( 𝑦) ,𝐿𝐿𝑚, 𝑛 = 1,2,3 ⋯(2.31)を与えた。また、式(2.11)を時間離散化して「反応項」およびゲージ場存在下での従属
変数𝛹𝑖,𝑗 (オーダーパラメータ)の更新式、
𝑖,𝑗

𝑖−1,𝑗

′
𝛹𝑖,𝑗
= 𝑎𝑖,𝑗 𝛹𝑖 + 𝑏𝑖,𝑗 (𝑤𝑥 𝛹𝑖+1,𝑗 + 𝑤
̅𝑥

𝑖,𝑗

𝑖,𝑗−1

𝛹𝑖−1,𝑗 + 𝑤𝑦 𝛹𝑖,𝑗+1 + 𝑤
̅𝑦𝛹𝑖,𝑗−1 )(2.32)′を得た。式(2.32)において𝛹𝑖,𝑗
は𝛹𝑖,𝑗 が時間𝜏後に更新された値である。

ただし、
𝑎𝑖,𝑗 = exp (−𝑏𝑖,𝑗 =

𝜎𝑖,𝑗 𝜏
)ℎ2𝛾(2.33)

1 − 𝑎𝑖,𝑗𝜎𝑖,𝑗(2.34)である。ここで、
𝒓𝑖,𝑗 + 𝒓𝑖+1,𝑗
𝒓𝑖,𝑗 + 𝒓𝑖,𝑗+1
𝑖,𝑗𝑖,𝑗𝑤𝑥 = exp (−iℎ𝐴𝑥 ()), 𝑤𝑦 = exp (−iℎ𝐴𝑦 ())
22である。ここで𝒓𝑖,𝑗 は格子点(𝑖, 𝑗)を示す位置ベクトルである。また、
̃𝑖,𝑗 |2 )ℎ2
𝜎𝑖,𝑗 = 4 + (𝛼 + 𝛽|𝛹

(2.35)(2.36)である。ここで̃𝑖,𝑗 =
𝛹

1 𝑖,𝑗
𝑖−1,𝑗
𝑖,𝑗
𝑖,𝑗−1
(𝑤 𝛹
+𝑤
̅𝑥 𝛹𝑖−1,𝑗 + 𝑤𝑦 𝛹𝑖,𝑗+1 + 𝑤
̅𝑦 𝛹𝑖,𝑗−1 )4 𝑥 𝑖+1,𝑗(2.37)である。
続いて、ゲージ場に関する方程式に関する離散化を以下の通り実装した。式(2.12)、(2.13)について、
𝑖,𝑗

𝑖,𝑗

𝑤𝑥 = exp(i𝜃𝑥 ) ,

𝑖,𝑗

𝑖,𝑗

𝑤𝑦 = exp (i𝜃𝑦 )(2.38)

d 𝑖,𝑗
𝜃
d𝑡 𝑦

(2.39)

および
𝑖,𝑗

𝜔𝑥 =

d 𝑖,𝑗
𝜃 ,
d𝑡 𝑥

𝑖,𝑗

𝜔𝑦 =

とおく。すると次の関係式が得られる。
𝑖,𝑗𝜔𝑦 = −𝑖,𝑗𝜔𝑦 = −

11
̅𝑖,𝑗 𝑤𝑥𝑖,𝑗 𝛹𝑖+1,𝑗 + 𝐶𝑧̅ 𝑖,𝑗−1 𝐶𝑧𝑖,𝑗 ]
Im[𝛹
𝜏𝐴
ℎ2

(2.40)

11
̅𝑖,𝑗 𝑤𝑦𝑖,𝑗 𝛹𝑖,𝑗+1 + 𝐶𝑧̅ 𝑖,𝑗 𝐶𝑧𝑖−1,𝑗 ]
Im[𝛹
𝜏𝐴
ℎ2

(2.41)

これらの関係式より、次のような、リンク変数の偏角に関するオイラー法が考えられる。

- 27 -

2.5

𝜏 𝑖,𝑗
𝑖,𝑗
𝑖,𝑗
𝑖,𝑗
𝜃𝑥 (𝑡 + 𝜏) = 𝜃𝑥 (𝑡) + (𝜔𝑥 (𝑡) + 𝜔𝑥 (𝑡 + 𝜏)) ,2(2.42)

𝜏 𝑖,𝑗
𝑖,𝑗
𝑖,𝑗
𝑖,𝑗
𝜃𝑦 (𝑡 + 𝜏) = 𝜃𝑦 (𝑡) + (𝜔𝑦 (𝑡) + 𝜔𝑦 (𝑡 + 𝜏))
2(2.43)

電界、磁界、電流密度の実装

TDGL 方程式(1.10)、(1.11)は次のようにも表される。
𝜕𝛾 ( + i𝑞𝜙) 𝛹 = (∇ − i𝑞𝑨)2 𝛹 − 𝛼𝛹 − 𝛽|𝛹 |2 𝛹𝜕𝑡

(2.44)∂𝑨̅ (∇ − i𝑞𝑨)𝛹] − ∇ × ∇ × 𝑨
γ ( + ∇𝜙) = 𝑞Im[𝛹∂𝑡(2.45)

ここで、𝜙はスカラーポテンシャル、𝑞は電荷である。
電界𝑬および磁界𝑩は次のように表される。
𝑬 = −∇𝜙 −

𝜕𝑨
,
𝜕𝑡

𝑩=∇×𝑨

(2.46)

ここで、式(2.44)、(2.45)、(2.46)はゲージ対称性を持つ。すなわち、任意関数𝜒のゲージ
変換
𝑨 → 𝑨 + ∇𝜒,

𝜙→𝜙−

𝜕𝜒
,
𝜕𝑡

𝛹 → 𝛹exp (i𝑞𝜒)(2.47)

によって式(2.44)、(2.45)、(2.46)は形を変えない。
このゲージ対称性を利用して、ゲージを選択する。今回はスカラーポテンシャル𝜙があらゆる場所で常にゼロであるゲージを採用した。このゲージのもとでは、式(2.44)、(2.45)、(2.46)はそれぞれ次のように表される。
∂𝛹= (∇ − i𝑞𝑨)2 𝛹 − 𝛼𝛹 − 𝛽|𝛹|2 𝛹∂𝑡(2.48)∂𝑨̅ (∇ − i𝑞𝑨)𝛹] − ∇ × ∇ × 𝑨= 𝑞Im[𝛹∂𝑡(2.49)

𝛾
𝜏𝐴

𝑬=−

∂𝑨
,
∂𝑡

𝑩= ∇×𝑨(2.50)

2.5.1 電界式(2.50)の左式および式(2.7)より、電界𝑬 = (𝐸𝑥 , 𝐸𝑦 )は次のように計算され、これらの式を実装した。

- 28 -

𝑖,𝑗
1 d𝜃𝑥
𝑖,𝑗
𝑖,𝑗
𝐸𝑥 =
= 𝜔𝑥 /ℎ(2.51)

𝑖,𝑗
1 d𝜃𝑦
𝑖,𝑗
𝑖,𝑗
𝐸𝑦 =
= 𝜔𝑦 /ℎ(2.52)ℎ d𝑡ℎ d𝑡

2.5.2 磁界𝑖,𝑗式(2.50)右式および式ループ変数𝐶z の定義式より、磁界は
𝑖,𝑗𝐵𝑧 = −

1𝑖,𝑗arg [𝐶𝑧 ]ℎ2(2.53)𝑖,𝑗と表され、これを実装した。磁界はループ変数𝐶z の値が定義されている空間位置と同じ位置で定義される。

2.5.3 電流密度式(2.49)右辺の第 1 項は電流密度に対応する。これを𝑱とおく。
̅ (∇ − i𝑞𝑨)𝛹]𝑱 = 𝑞Im[𝛹(2.54)であり、𝑖,𝑗
̅𝑖,𝑗 𝑤𝑥𝑖,𝑗 𝛹𝑖+1,𝑗 ]/ℎ𝐽𝑥 = 𝑞Im[𝛹𝑖,𝑗̅𝑖,𝑗 𝑤𝑦𝑖,𝑗 𝛹𝑖,𝑗+1 ]/ℎ
𝐽𝑦 = 𝑞Im[𝛹(2.55)(2.56)

である。式(2.55)、式(2.56)を実装した。電流密度の𝑥成分、𝑦成分は対応するリンク変数𝑤𝑥 および𝑤𝑦 が定義されている空間位置と同じ位置で定義される。

- 29 -

2.6

ピンの導入

ピンの導入に関しては、凝縮エネルギーの差により生じる引力相互作用を想定して実装した。具体的には TDGL 方程式中のパラメータ𝛼に空間依存性を持たせることによって実現した。超伝導領域において𝛼 = −20としているため、ピンを配置したい座標で𝛼 >−20とすることで超伝導部分とその位置との間に凝縮エネルギー的な差が生じるためその部分は磁束線を固定するピンとして働く。

2.7

描画

描画処理は統合開発環境 Processing により行った。10 𝜉 × 10 𝜉の大きさを持つ超伝導
平面について、左右 2 つの画面を用意した。左の画面にはオーダーパラメータに関する
画面を割り当てた。各部分のオーダーパラメータの大きさの2乗(|𝛹 |2 )を輝度に、位相を色相に割り当て表示させた。これにより量子化磁束の観測を行った。ピンを配置した領域には白い格子を描画し、ピンの位置が目に見えるように実装した。また、右の画面には、磁界の大きさを色相で、電界を黒矢印で、電流密度を白矢印で表示させた。これにより電磁現象の可視化を行った。

2.8

リザバーコンピューティング

リザバーコンピューティングのタスクを行う前段階として、描画された超伝導領域か
らの電界の計算を行った。Processing 上で描画された大きさ10 𝜉 × 10 𝜉２次元超伝導領域に電流密度の時間変化を与え、描画領域から無作為に50個の点を選択し、選択した点
における電界を時間刻み幅𝜏=0.001 ごとに計算した。ただし、シミュレーション開始からしばらくは過渡状態であり安定した入出力応答が得られないため、定常状態に入ってから電界の計算を行った。これを𝑡 = 𝑡1 ~𝑡10000 まで行った。
計算した電界値𝐸(𝑡𝑗 )を入力信号系列である電流密度𝐽(𝑡𝑗 )に駆動されるリザバーのノード𝑿(𝑡𝑗 )とした。この電流密度と電界を用いてリザバーコンピューティングのタスクくを行った。Fig. 2.5 に電流密度の印加と電界の抽出に関する概略図を示す。

- 30 -

Fig. 2.5 電流の印加方法と電界の抽出方法の概略図
続いて実際に実行したタスクについて詳細に説明する。卒業研究ではリサージュ波形と FFT による入出力応答の調査、波形生成タスク、NARMA2 タスク、非線形-メモリタ
スクの主に 4 つのタスクをそれぞれ 1 回ずつ行った。入出力応答の調査では超伝導がリザバーが持つべき性質を備えていることが明らかにされた。また、波形生成タスクではいずれの波形に対しても 8 割以上高精度での予測が実現された。本研究は精度向上を目的としているため、タスクは精度の向上が見られた場合わかりやすい NARMA2 タスクと非線形-メモリタスクの 2 つに絞った。ここでそれぞれのタスクの詳細を説明する。
NARMA2 タスクでは、全 10000 データのうち𝑡500 ~𝑡800 までの 300 データを利用した。
この 300 データを8 ∶ 2に分割して、前半 240 データを学習、後半 60 データを予測に用いた。教師信号として式(1.91)を用いた。この式において𝑚 = 2、𝑎1 = 0.4、𝑎2 = 0.4、
𝑎3 = 0.6、𝑎4 = 0.1とした。学習には過学習を抑制するためにリッジ回帰(式(1.90))を用いた。性能の評価に関しては、評価指標として R2(式(1.92))を用いた。
非線形メモリタスクでは、全 10000 データのうち𝑡5000 ~𝑡5500 の 500 データを利用した。学習部分と予測部分の割合は NARMA2 タスク同様8 ∶ 2とした。教師信号として式(1.94)を用いた。評価指標には R2 を用いた。式(1.94)の𝜈と𝜏について𝜈は 0 から 0.1 ずつ
5.0 まで、𝜏は 0 から 1 ずつ 10 まで値を変化させ、それぞれの場合についての R2 を計算した。
これら 2 つのタスクを超伝導リザバーのパラメータを変化させながら行った。以下にその詳細を記述する。

- 31 -

2.8.1 ピンの数の変更
まず、超伝導リザバーに入れるピンの数を変化させた。具体的には、ピンの数を 0 から 50 まで 10 ずつ変化させた。このときピンを入れる位置はランダムとし、ピン 1 つあたりのピン力の強さに関係するパラメータ𝛼は 0 に固定した。た。それぞれのピンの数において NARMA2 タスクを行い、その結果を横軸に時間、縦軸に出力の値をとり、図に目標出力、学習時の出力、予測時の出力に関して色分けして表した。また、ピンの数
を横軸、R2 を縦軸にとり、ピンの数の変化に対する R2 の変化に関しても図に示した。
また、非線形-メモリタスクを行いそれぞれの場合において、𝜈を固定した状態での𝜏の変化に対する R2 の変化、𝜏を固定した状態での𝜈の変化に対する R2 の変化を図に示した。

2.8.2 ピンの配置の変更
続いて、ピンの配置に変化を与えた。ピンを左右の境界付近、上下の境界付近、中心付近に配置した場合の 3 つのパターンを用意しそれぞれに関して NARMA2 タスクと非線形-メモリタスクを行った。このとき、ピンの数は 10 に固定し、ピン力の強さに関係するパラメータ𝛼は 0 に固定した。NARMA2 タスクに関して、ピンの数を変化させた場
合同様に、横軸に時間、縦軸に出力の値をとり、図に目標出力、学習時の出力、予測時の出力に関して色分けして表した。また、それぞれの場合の R2 を表にまとめた。非線形-メモリタスクに関して、ピンの数を変化させた場合同様に結果を図に示した。

2.8.3 ピン 1 つあたりのピン力の変更
最後に、ピン 1 つあたりのピン力を変化させるためにパラメータ𝛼を変化させた。超伝導リザバー全体としては𝛼 = −20に固定したままで、ピンがある位置の𝛼を 0 から 20まで 5 ずつ変化させた。このとき、ピンの数は 30 に固定し、配置はランダムとした。
𝛼の値それぞれに関して NARMA2 タスクと非線形-メモリタスクを行った。結果はピンの数を変化させた場合同様に図に示した。

- 32 -

第3章
3.1

結果と考察

量子化磁束および電磁現象の可視化

各条件の変化に対する量子化磁束の運動およびそれに伴う電磁現象の可視化の結果を以下に示す。

3.1.1 ピンの数を変化させた場合
ピンの数を変化させていった場合での磁束量子の動きと電磁現象を可視化したもの
を Fig. 3.1 に示す。本来は動画であるが、𝑡 = 𝑡5000 のときを画像として切り取っている。
まず、ピンの数に関わらず共通して確認されたことについて述べる。左図が、オーダーパラメータの大きさと位相(磁束線の動き)に関して示した図である。オーダーパラメータの大きさを輝度、位相を色相で表している。白いメッシュはピンを表している。右
図が、電磁現象に関して示した図である。磁束密度を色相、電流密度を白矢印、電界を黒矢印で表している。
左図において、黒い部分で|𝛹 |2 = 0となっており、その部分に磁束線が侵入している。
侵入した磁束線はだいたい三角格子を構成していることが見て取れた。また、磁束線 1本に注目したとき、その周りを一周する円を考えたとき、赤を基準とすると赤から 1 周
してまた赤に戻る様子が見られた。つまり、位相が2πn(ただし、𝑛は自然数)ずれている
様子が確認された。シミュレーション開始から過渡状態の間、この磁束線が左端から侵入し右側(𝑥軸正の向き)に移動していく様子を観測することができた。これは、電流密度を𝑦軸正の向き、磁束密度を𝑧軸正の向きに与えたことで磁束線にローレンツ力𝐹𝑙 =𝐽 × 𝐵が働いたためであると考えられる。定常状態に入った後、与える電流密度は−1 ≤𝐽𝑦 ≤ 1の乱数で時間変化させたため、磁束線は𝐽𝑦 > 0のときに𝑥軸正の向き、𝐽𝑦 < 0のときに𝑥軸負の向きに運動した。また、磁束線の運動速度は|𝐽𝑦 |が大きいときには早く、小さいときには遅くなった。
右図において、磁束密度は磁束線の侵入している部分で大きくなっている様子が見られた。また、電流密度は全体で見ると𝑦軸正の向きに流れているが、磁束線侵入部分ではその周りを回るように流れている様子が確認された。磁束線の周りを回る電流が量子化磁束を構成する要素であると考えられる。また、電界に関しては𝐸 = 𝐵 × 𝑣にしたがって磁束線の運動と垂直に発生する様子を確認した。
続いてピンの数の変化によって見られた違いを説明する。ピンの数が 0 の場合には、
磁束線は止まることなくゆっくりと移動していく様子が確認された。ピンを入れた場合にはピンのある位置で磁束線がしばらくの間固定される様子が見られた。その結果磁束線が構成する三角格子がやや乱れるようになった。また、ピンの数が多くなるほど固定- 33 -

される磁束線の数も多くなり、ピンに止められていない磁束線の運動も鈍くなる様子が確認された。これは磁束線同士にばねのような弾性力働き磁束バンドル(1.5.2 小節参照)
を構成しているためであると考えられる。

- 34 -

Fig. 3.1 ピンの数ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。上から順に、ピンの数が 0、10、20、30、40、50 となっている。

3.1.2 ピンの配置を変化させた場合
ピンの配置を変化させた場合の量子化磁束の運動とそれに伴う電磁現象の可視化結果
を Fig.3.2 に示す。まず、中央付近に配置した場合には、磁束線はピンのある位置では- 35 -

ピンの並びに合わせて整列された状態になった。一度ピンに固定された磁束線は電流密度の変化に対して、ピンの中でわずかに左右に動くだけであった。ピンのある部分
を境にその左右の領域ではピンの近くでは磁束線の動きが鈍く、ピンから離れたとこ
ろではスムーズに動く様子が確認された。次に左右境界付近に配置した場合には、磁束線が侵入したときと磁束線が出ていくときに磁束線がピンに止められて、その次に来る磁束線に押し出される形でピンに止められている磁束線がピンから外れ内部に侵入、または外部に出る様子が見られた。最後に上下境界付近に配置した場合には、上下のピンに止められている磁束線以外はスムーズに運動し、止められている磁束線は
基本的には止められたままでたまにピンから外れる様子が確認された。また、通常磁
束線は左右境界にに出入りするが、上下境界の左端や右端から出入りする磁束線がい
くつか見られた。

- 36 -

Fig. 3.2 ピンの配置ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。上から順に、ピンを中央付近に配置した場合、左右境界付近に配置した場合、上下境界付近に配置した場合を表している。

3.1.3 ピン 1 つあたりのピン力を変化させた場合
ピン１つあたりのピン力に関係するパラメータ𝛼を変化させた場合の量子化磁束運動およびそれに伴う電磁現象の可視化結果を Fig. 3.3 に示す。𝛼 = 0のときには上までの結果と同様ピンは白いメッシュとして観察されたが、𝛼 = 5以降はピンが黒い正方形として観察された。また、𝛼が大きくなるにつれてピンはより濃い黒色になり𝛼 = 20で最も濃くなった。これは、本研究で想定されているピンが凝縮エネルギー相互作用を前提に実装されていることが原因である。このシミュレーションではオーダーパラメータの大きさ|𝛹 |2 が色相で表されている。|𝛹 |2 は𝛼に逆比例するため、𝛼を大きくすることで|𝛹 |2が小さくなり黒く表されるようになったと考えられる。
また、𝛼を大きくしていくと、磁束線がより長い間止められる様子が観察された。これは𝛼の増加によりピン 1 つあたりのピン力が大きくなり、磁束線に働くローレンツ力に対してピン力が釣り合っていられる時間が長くなったためであると考えられる。

- 37 -

- 38 -

Fig. 3.3 ピン 1 つあたりのピン力に関係するパラメータ𝛼ごとの量子化磁束の運動およびそれに伴う電磁現象の可視化。上から順に、𝛼 = 0, 5, 10, 15, 20のときを表して
いる。

3.2

NARMA2 タスク

各条件の変化に関する NARMA2 タスクの結果を以下に示す。

3.2.1 ピンの数を変化させた場合
ピンの数を変化させたときにおけるピンの数ごとの NARMA2 タスクの結果を Fig.
3.4 に示す。赤い線が目標出力、青い線が学習時の出力、緑の線が予測時の出力である。

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

y

y

いずれの場合においても出力は大まかに目標出力を再現している。

0

−0.2

−0.2

−0.4

−0.4

−0.6

500

−0.6

(a)

−0.8
550

600

650

700

750

(b)

−0.8
800

500

550

600

650

t

t

- 39 -

700

750

800

1
0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

y

y

1
0.8

0

−0.2

−0.2

−0.4

−0.4

−0.6

500

−0.6

(c)

−0.8
550

600

650

700

750

(d)

−0.8

800

500

550

600

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

−0.2

−0.4

−0.4

−0.6

−0.6

(e)

−0.8
550

600

650

700

750

800

0

−0.2

500

650

t

y

y

t

700

750

(f)

−0.8

800

t

500

550

600

650

700

750

800

t

Fig. 3.4 各ピンの数における NARMA2 タスクの結果。横軸が時間、縦軸が出力を表している。(a)はピンの数が 0、(b)はピンの数が 10、(c)はピンの数が 20、(d)はピンの数が 30、(e)はピンの数が 40、(f)はピンの数が 50 のときの結果を表している。

次に、ピンの数の変化に対する NARMA2 タスクの R2 の変化を Fig. 3.5 に示す。ピンの数が𝑥軸、R2 が𝑦軸に対応している。ピンの数が 0 から 10 に増えたときは R2 が増加したが、それ以降において R2 は減少した。これは、ピンの増加による超伝導リザバーの非線形性の空間的な差の変化が原因であると考えられる。ピンの数が 0 とはピンが無いということであるため、ピンの数が 0 から 10 に増えるとき、その変化はピン無の状
態からピン有の状態への変化と考えることができる。ピンが無い状態では超伝導リザバ
ー内の非線形性は空間的にほぼ一様である。その状態からピン有の状態に変化すること
でピンのある部分はピンの無い部分より非線形性が大きくなり、空間的な非線形性の差が生じる。そのため R2 が増加したと考察する。また、ピンの数を 10 から増やすとき、
動く磁束線とピンに止められる磁束線が理想的な塩梅で分かれていた状態から、磁束線が全体的に止められる状態となり、非線形性の空間的な差が小さくなる。その結果 R2は減少したと考えられる。

- 40 -

1

0.8

R

2

0.6

0.4

0.2

0
0

10

20

30

40

50

Number of pins
Fig. 3.5 電流密度の乱数波入力に対する電界の出力応答

3.2.2 ピンの配置を変化させた場合
ピンの配置を変化させたときに関する各ピンの配置における NARMA2 タスクの結果を Fig. 3.6 に示す。軸および色の対応は Fig. 3.4 と同様である。こちらも目標出力の時系列変化をおおまかに予測できている様子を確認できた。
1
0.8
0.6
0.4

y

0.2
0
−0.2
−0.4
−0.6

(a)

−0.8
500

550

600

650

700

750

800

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

y

y

t

0

−0.2

−0.2

−0.4

−0.4

−0.6

500

−0.6

(b)

−0.8
550

600

650

700

750

(c)

−0.8

800

500

t

550

600

650

700

750

800

t

Fig. 3.6 各ピンの配置における NARMA2 タスクの結果。(a)がピンを中央付近に配置
した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付近に配置した場合の結果である。

- 41 -

次に各ピンの配置に対する R2 について Table 3.1 に示す。1 列目にピンの配置、2 列
目に R2 を示している。いずれのピンの配置においても R2 はおおよそ 0.3 程度となった。ピンの配置に対して R2 の変化がほとんど生じなかった原因として、ピンの数の空間的偏りに変化が無かったことが考えられる。まず、ピンの数はいずれの場合も 10 に固定されている。また、ピン 1 つあたりのピン力もすべてのピンにおいて等しい。今
回のピンの配置はいずれも境界付近や中央付近といった極端な場所に集中している。
したがって超伝導リザバー内は、いずれの場合も磁束線がまったく止められない領域とピンが集中していて磁束線がよく止められる領域に 2 極化するという共通した状態になっている。そのため非線形の空間的な差が同程度となり R2 の差も小さくなったと考察する。
Table. 3.1 ピンの配置に対する R2
Pin Placement

Accuracy R2Center0.319

Up and Down

0.327

Left and Right

0.301

3.2.3 ピン 1 つあたりのピン力を変化させた場合
ピン 1 つ当たりのピン力に関係するパラメータ𝛼を変化させた場合の NARMA2 タス
クの結果を Fig. 3.7 に示す。軸及び線の色の対応は Fig. 3.4、Fig. 3.6 と同様である。この場合においてもおおよその目標出力へ対しての予測はできており、見た目においてはパラメータ𝛼の違いによる結果の違いは見られなかった。
1
0.8
0.6
0.4

y

0.2
0
−0.2
−0.4
−0.6

(a)

−0.8
500

550

600

650

t

- 42 -

700

750

800

1
0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

y

y

1
0.8

0

−0.2

−0.2

−0.4

−0.4

−0.6

500

−0.6

(b)

−0.8
550

600

650

700

750

(c)

−0.8
500

800

550

600

1
0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

y

y

1
0.8

700

750

800

0

−0.2

−0.2

−0.4

−0.4
−0.6

−0.6

(d)

−0.8
500

650

t

t

550

600

650

700

750

(e)

−0.8
800

500

550

600

650

700

750

800

t

t

Fig. 3.7 各𝛼の値に対する NARMA2 タスクの結果。(a)が𝛼 = 0、(b)が𝛼 = 5、(c)が
𝛼 = 10、(d)が𝛼 = 15、(e)が𝛼 = 20に対応している。

続いて、パラメータ𝛼の変化に対する NARMA2 タスクの R2 の変化を Fig. 3.8 に示す。𝛼 = 0から𝛼 = 5としたとき R2 は減少し、それ以降において R2 は増加する結果になった。このような結果となった原因として、まず𝛼 = 0とそれ以外について考える。
𝛼 = 0のピンは、ピンの数を変えた場合、ピンの配置を変えた場合のときに入れたものと全く同じものである。可視化シミュレーションでは白いメッシュとして表されていた。𝛼 = 0のピンはちょうど超伝導状態と常伝導状態の境の状態にあると言える。一方
それ以外のピンは可視化シミュレーションでは黒い正方形として表されていたが、こちらは完全に常伝導状態にあると言える。つまり、𝛼 = 0のピンはそれ以外のピンと比較して不安定なピンであるため複雑な非線形性を生み、高めの R2 を得ることができたのではないかと考える。また、それ以外のピンにおいて𝛼が大きくなるにつれて R2 が大きくなったが、これに関しては𝛼の増加とともにピンある部分と無い部分の凝縮エネ
ルギーの差が大きくなってより強い非線形性を得ることができ R2 の上昇につながったと考えられる。

- 43 -

1

0.8

R

2

0.6

0.4

0.2

0
0

5

10



15

20

Fig. 3.8 入出力データの FFT

3.3

非線形-メモリタスク

各条件の変化に関する非線形-メモリタスクの結果を以下に示す。

3.3.1 ピンの数を変化させた場合
非線形-メモリタスクのピンの数を変化させたときの各結果において、メモリパラメータ𝜏 = 1と固定して非線形性パラメータ𝜈を変化させた場合の結果を Fig. 3.9 に示す。
𝑥軸が非線形性パラメータ𝜈、𝑦軸が R2 を表している。いずれの場合も𝜈の値が 2 辺りから R2 が減少し始め𝜈 = 4から𝜈 = 4.5の間辺りで R2 がほぼゼロになる様子が確認された。
目視では違いが分かりにくいため、𝜈 = 3.0での R2 の値を比較した。ピン 0 個ではR2 =0.72、ピン 10 個ではR2 = 0.74、ピン 20 個ではR2 = 0.74、ピン 30 個ではR2 = 0.74、ピン 40 個ではR2 = 0.65、ピン 50 個ではR2 = 0.69となった。ピンの個数が 0 から 10 に増えるとき R2 が増加し、それ以上ピンの数を増やすと R2 が減少するという NARMA2 タ
スクの R2 のピンの個数依存性の結果におおむね従う結果となった。したがって、
NARMA2 タスクのピンの数の変化に対する R2 の変化は純粋に非線形性の変化によるものであったと考えられる。

- 44 -

0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(a)

(b)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5





0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(d)

(c)
0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5





0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(e)

(f)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
2

Fig. 3.9 各ピンの数における非線形-メモリタスクの R の𝜈依存性。(a)がピン 0 個、(b)がピン 10 個、(c)がピン 20 個、(d)がピン 30 個、(e)がピン 40 個、(f)がピン 50 個のときの結果を表している。

- 45 -

続いて、非線形-メモリタスクのピンの数を変化させたときの各結果において、非線
形性パラメータを𝜈 = 0.1と固定してメモリパラメータ𝜏を変化させた場合の結果を Fig.
3.10 に示す。𝑥軸がメモリパラメータ𝜏、𝑦軸が R2 を表している。いずれの場合も𝜏の値の増加とともに R2 が減少する様子が確認された。ピンの数が 0 のときは、それ以外のときに比べて𝜏 = 2の R2 が低くなった。また、ピンの数が 50 のときの結果はそれ以外の結果よりも全体的に記憶特性が高くなった。ピンの数が 0 の結果からピンの数が 50の結果まで順番に見ていくと、やや上下はあるもののピンの数の増加とともに記憶特性も増加しているように見て取れる。したがって、ピンの数は超伝導リザバーの記憶特性に影響を与えるということが明らかになった。

0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0
1

0.2

(a)
2

(b)
3

4

5

6

7

8

0
1

9

2

3

4



5

6

7

8

9

6

7

8

9



0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(d)

(c)
0
1

2

3

4

5

6

7

8

9



0
1

2

3

4

5



- 46 -

0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(e)
0
1

2

(f)
3

4

5

6

7

8

0
1

9

2

3



4

5

6

7

8

9


2

Fig. 3.10 各ピンの数における非線形-メモリタスクの R の𝜏依存性。(a)がピン 0 個、(b)がピン 10 個、(c)がピン 20 個、(d)がピン 30 個、(e)がピン 40 個、(f)がピン 50 個のときの結果を表している。

3.3.2 ピンの配置を変化させた場合
非線形-メモリタスクの結果において、メモリパラメータ𝜏 = 1と固定して非線形性パ
ラメータ𝜈を変化させた場合のピンの配置ごとの結果を Fig. 3.11 に示す。ピンの数を変化させた場合と同様にいずれの場合も𝜈の値が 2 辺りから R2 が減少し始め𝜈 = 4から𝜈 =4.5の間辺りで R2 がほぼゼロになる様子が確認された。こちらの場合も𝜈 = 3.0のときのR2 を比較してみると、ピンを中央付近に配置したときR2 = 0.62、左右境界に配置した
ときR2 = 0.66、上下境界付近に配置したときR2 = 0.66であった。つまり R2 はほとんど変化していないと言える。したがって、ピンの配置によって非線形性に大きな変化は生じないことがわかった。これは NARMA2 タスクの結果とも一致している。
1

0.8

R

2

0.6

0.4

0.2

(a)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

- 47 -

0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0.2

(b)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5(c)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
2

Fig. 3.11 各ピンの配置における非線形-メモリタスクの R の𝜈依存性。(a)がピンを中央付近に配置した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付近に配置した場合の結果を表している。

続いて、非線形-メモリタスクのピンの配置を変化させたときの各結果において、非線形性パラメータを𝜈 = 0.1と固定してメモリパラメータ𝜏を変化させた場合の結果をFig. 3.12 に示す。いずれの結果も𝜏の増加とともに R2 の減少傾向が見られた。しかしながら、3 つの結果には現象の仕方に違いが見られた。ピンを中央付近に配置した場合においては、𝜏 = 4で R2 が 0.3 程度に落ちた後ほぼ横ばいで微減していく結果となった。
ピンを左右境界に配置した場合においては、R2 はピンを中央付近に配置した場合に比べて緩やかに減少し、𝜏 = 9で R2 が約 0.2 付近に落ち着いた。ピンを上下境界付近に配置した場合においては、𝜏 = 1から𝜏 = 3まではほかの場合と比べて高い R2 を保っていたがが、𝜏 = 4以降も減少を続け𝜏 = 9では R2 はほぼ 0 まで落ちた。まず、ピンの配置が上下境界付近のときのみ R2 がほぼ 0 まで減少したことについて、ピンの整列の方向が原因であると考える。ピンの整列の方向は上下境界付近のときのみ横向き、それ以外のときは縦向きである。磁束線の運動の方向は横向きであるため、より良く磁束線を止めることができるのは縦方向に整列したピンである。整列が横方向だとピンに止められない
磁束線はほぼピンの影響を受けることなく運動できるため、記憶特性に関与できるのは精々ピンの周辺の領域のみである。反対に整列が縦方向の場合、ピンに止められない磁束線も進行方向または逆方向のすぐ近くにピンに止められた磁束線があると、その磁束線から影響を受ける。また、その近くの近くの磁束線も影響を受け超伝導リザバー全体
が記憶特性に関与することができる。そのためピンを中央付近に配置した場合と左右境界に配置した場合で𝜏が大きくなってもある程度の R2 の大きさを保っていたのだと考えられる。また、ピンを左右境界に配置した場合の方がピンを中央付近に配置した場合- 48 -

よりも全体的に R2 の減少が緩やかになったことについては、場所による磁束線の運動の激しさが関係している考える。磁束線の運動は中央付近よりも左右境界付近の方が激しい。したがって、磁束線がピンに止められることによる影響はピンを左右境界に配置したときの方が大きくなる。そのためピンを左右境界に配置した場合の方が全体的な記憶特性が良くなったのだと考えられる。

1

0.8

R

2

0.6

0.4

0.2

(a)

0
1

2

3

4

5

6

7

8

9



0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0
1

0.2

(b)
2

3

4

5

6

7

8

9

0
1

(c)
2

3



4

5

6

7

8

9



Fig. 3.12 各ピンの配置における非線形-メモリタスクの R2 の𝜏依存性。(a)がピンを中央付近に配置した場合、(b)が左右境界付近に配置した場合、(c)が上下境界付近に配置した場合の結果を表している。

- 49 -

3.3.3 ピン 1 つあたりのピン力を変化させた場合
非線形-メモリタスクの結果において、メモリパラメータ𝜏 = 1と固定して非線形性
パラメータ𝜈を変化させた場合のパラメータ𝛼ごとの結果を Fig. 3.13 に示す。上 2 つの場合の結果同様いずれの場合も𝜈の値が 2 辺りから R2 が減少し始め𝜈 = 4から𝜈 = 4.5の間辺りで R2 がほぼゼロになる様子が確認された。𝜈 = 3.0のときのそれぞれの R2 は、
𝛼 = 0でR2 = 0.65、𝛼 = 5でR2 = 0.68、𝛼 = 10でR2 = 0.70、𝛼 = 15でR2 = 0.70、𝛼 = 20
でR2 = 0.69となった。𝛼 = 0のときの R2 のみやや低い結果になった。純粋な非線形性のみを見ると NARMA2 タスクで R2 が大きかった𝛼 = 0のときよりも𝛼 = 5や𝛼 = 10の方が非線形性が強いと言える。したがって NARMA2 タスクの結果は純粋な非線形性以外の要素が関わっていると考えられる。

1

0.8

R

2

0.6

0.4

0.2

(a)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5



0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.2

0.4

0.2

(b)

(c)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5- 50 -

0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.2

0.4

0.2

(d)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5(e)

0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
2

Fig. 3.13 各𝛼の値における非線形-メモリタスクの R の𝜈依存性。(a)が𝛼 = 0の場合、(b)が𝛼 = 5の場合、(c)が𝛼 = 10の場合、(d)が𝛼 = 15の場合、(e)が𝛼 = 20の場合の結果を表している。

続いて、非線形-メモリタスクのピン 1 つあたりのピン力に関係するパラメータ𝛼を変化させたときの各結果において、非線形性パラメータを𝜈 = 0.1と固定してメモリパ
ラメータ𝜏を変化させた場合の結果を Fig. 3.14 に示す。他の条件変更の場合と同様、𝜏の増加に対して R2 は減少傾向となった。𝛼 = 0のときの R2 の減少傾向が最も緩やかであり全体的に見た記憶特性も最も良い結果となった。それ以外の𝛼の場合の結果を比較しても、𝛼の値の増加に対しての記憶特性の向上または低下は見られなかった。このことを𝜈依存性を調べたときの結果と比較して考えると、非線形性が比較的弱い𝛼 = 0のとき記憶特性が向上し、非線形性が比較的強いそれ以外の𝛼のとき記憶特性が低下したと捉えることができる。一般にリザバーコンピューティングにおいて非線形性と記憶性はトレードオフの関係にあると言われており、この結果はそれに従っていると言える。𝛼の違いごとに生じているわずかな記憶特性の差異は、ピンをランダムに配置したことによる違いであると考える。

- 51 -

1

0.8

R

2

0.6

0.4

0.2

(a)

0
1

2

3

4

5

6

7

8

9



0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0
1

0.2

(b)
2

3

4

5

6

7

8

(c)

0
1

9

2

3

4



5

6

7

8

9

6

7

8

9



0.8

0.8

0.6

0.6

R

R

2

1

2

1

0.4

0.4

0.2

0
1

0.2

(d)
2

3

4

5

6

7

8

9

(e)

0
1

2



3

4

5


2

Fig. 3.14 各𝛼の値における非線形-メモリタスクの R の𝜏依存性。(a)が𝛼 = 0の場合、(b)が𝛼 = 5の場合、(c)が𝛼 = 10の場合、(d)が𝛼 = 15の場合、(e)が𝛼 = 20の場合の結果を表している。

- 52 -

第4章

まとめ

本研究では、リザバーコンピューティングにおいてリザバーの内部状態と時系列予測の結果との関係が不明瞭であるという課題に関して、その課題へのアプローチの１つとして、第 2 種超伝導体をもとにした超伝導リザバーを用いてその内部状態を可視化した。また、内部状態を変化させながら時系列予測のタスクを行いその内部状態の変化が予測精度に与える影響について調査した。
まず、Processing というソフトを用いて TDGL 方程式を AFI 法で数値的に解き 2 次元
領域で超伝導リザバーの可視化を行った。その際ピンの数、ピンの配置、ピン 1 つあたりのピン力を変化させながらそれぞれの場合においてシミュレーションした。超伝導体内での磁束線の振る舞い、およびそれに伴う超伝導体内の磁束密度、電流密度、電界などの様子を観察することができた。また、条件の違いによる磁束線の運動の変化、生じる電磁現象の変化を観察することができた。
続いてリザーバーコンピューティングのタスクの 1 つである NARMA2 タスクを、ピ
ンの数、ピンの配置、ピン 1 つあたりのピン力をそれぞれ変化させた場合について行った。ピンの数を変化させた場合には、ピンの数を 0 から 10 に増やしたとき予測精度 R2は減少し、それ以上ピンの数を増やすと R2 は減少していくという結果を得た。また、ピンの配置を変化させた場合には、ピンの配置によって予測精度はほとんど変化しないという結果を得た。また、ピン 1 つあたりのピン力に関係するパラメータ𝛼を変化させた場合には、𝛼 = 0から𝛼 = 5に変化させたとき予測精度は減少し、それ以降において予測精度が上昇するという結果を得た。
最後に直接リザバーの非線形性と記憶特性を測るタスクである非線形-メモリタスク
を NARMA2 タスクのとき同様、ピンの数、ピンの配置、ピン 1 つあたりのピン力をそれぞれ変化させた場合について行った。ピンの数を変化させた場合、非線形性について
はピンの数を 0 から 10 に増やしたときに非線形性が強くなり、さらにピンの数を増やすと非線形性が弱くなっていくという結果を得た。記憶特性に関しては、ピンの数を増やすと記憶特性も向上するという結果を得た。ピンの配置を変化させた場合、非線形性には大きな変化は見られなかった。記憶特性に関しては、ピンを左右境界に配置したときが最も高く、次いで中央付近に配置したとき、最も低かったのが上下境界付近に配置したときであった。ピン 1 つあたりのピン力を変化させた場合、非線形性については
𝛼 = 0のときのみ低くなり、他の𝛼の値では同程度であった。記憶特性に関しては、𝛼 =0のときのみ高くなり、それ以外で同程度であった。
本研究を通して超伝導リザバーの精度向上の方法を明らかにすることができた。今後は他の物理リザバーとも比較し、さらなる精度向上へ向けたアプローチ方法を探り、リザバーコンピューティング全体の普遍的な精度向上へと繋げていくことが重要である。

- 53 -

参考文献[1] Severino Segato, Giorgio Marchesini, Luisa Magrin, Barbara Contiero, Igino Andrighetto, andLorenzo Serva, A machine learning-based assessment of maize silage dry matter losses by
net-bags buried in farm bunker silos, Agriculture, 12 (2022) 785.
[2] Asif Mahmood, Ahmad Irfan, and Jin-Liang Wang, Machine learning and molecular dynamics
simulation-assisted evolutionary design and discovery pipeline to screen efficient small
molecule acceptors for ptb7-th-based organic solar cells with over 15% efficiency, Journal of
Materials Chemistry A, 10 (2022) 4170–4180.
[3] George Em Karniadakis, Ioannis G Kevrekidis, Lu Lu, Paris Perdikaris, Sifan Wang, and Liu
Yang, Physics-informed machine learning, Nature Reviews Physics 3 (2021) 422–440.
[4] Asmaa Abbas, Mohammed M Abdelsamea, and Mohamed Medhat Gaber, Classification ofcovid-19 in chest x-ray images using detrac deep convolutional neural network, AppliedIntelligence, 51 (2021) 854–864.
[5] Kohei Nakajima, Physical reservoir computing—an introductory perspective, Jpn. J. Appl.
Phys., 59 (2020) 060501.
[6] Herbert Jaeger, Echo state network, scholarpedia 2 (2007) 2330.
[7] Ali Rodan and Peter Tino, Minimum complexity echo state network, IEEE transactions onneural networks, 22 (2010) 131–144. 18[8] Herbert Jaeger, Tutorial on training recurrent neural networks, covering bppt, rtrl, ekf andthe ”echo state network” approach, GMD-Forschungszentrum Informationstechnik, 195(2002) 1–48.
[9] Mingzhe Chen, Walid Saad, and Changchuan Yin, Liquid state machine learning for resource
and cache management in lte-u unmanned aerial vehicle (uav) networks, IEEE Transactions
on Wireless Communications 18 (2019) 1504–1517.
[10] Mingzhe Chen, Walid Saad, and Changchuan Yin, Liquid state machine learning for resource
allocation in a network of cache-enabled lte-u uavs, In GLOBECOM 2017-2017 IEEE
Global Communications Conference, (2017) 1–6.
[11] Gouhei Tanaka, Toshiyuki Yamane, Jean Benoit H´eroux, Ryosho Nakane, Naoki
Kanazawa, Seiji Takeda, Hidetoshi Numata, Daiju Nakano, and Akira Hirose, Recent
advances in physical reservoir computing: A review, Neural Networks, 115 (2019) 100–
123.
[12] Yuki Usami, Bram van de Ven, Dilu G Mathew, Tao Chen, Takumi Kotooka, Yuya
Kawashima, Yuichiro Tanaka, Yoichi Otsuka, Hiroshi Ohoyama, Hakaru Tamukoh, et al,- 54 -

In-materio reservoir computing in a sulfonated polyaniline network, Advanced Materials,
33 (2021) 2102688.
[13] Wencong Jiang, Lina Chen, Kaiyuan Zhou, Liyuan Li, Qingwei Fu, Youwei Du, and RH
Liu, Physical reservoir computing using magnetic skyrmion memristor and spin torquenano-oscillator, Applied Physics Letters, 115 (2019) 192403.
[14] Yanan Zhong, Jianshi Tang, Xinyi Li, Bin Gao, He Qian, and Huaqiang Wu, Dynamic
memristor-based reservoir computing for high-efficiency temporal signal processing,Nature communications, 12 (2021) 1–9.
[15] Sumito Tsunegi, Tomohiro Taniguchi, Kohei Nakajima, Shinji Miwa, Kay Yakushiji, Akio
Fukushima, Shinji Yuasa, Hitoshi Kubota, Physical reservoir computing based on spin
torque oscillator with forced synchronization, Appl. Phys. Lett., 114 (2019) 164101.
[16] Yuichiro Yada, Shusaku Yasuda, Hirokazu Takahashi, Physical reservoir computing with
FORCE learning in a living neuronal culture, Appl. Phys. Lett., 119 (2021) 173701.
[17] Dan A. Allwood, Matthew O. A. Ellis, David Griffin, Thomas J. Hayward, Luca
Manneschi, Mohammad F. KH. Musameh, Simon O’Keefe, Susan Stepney, Charles
Swindells, Martin A. Trefzer, Eleni Vasilaki, Guru Venkat, Ian Vidamour, Chester Wringe,
A perspective on physical reservoir computing with nanomagnetic devices, Appl. Phys.
Lett., 122 (2023) 040501[18] Gerasimos Angelatos, Saeed A. Khan, and Hakan E. T¨ureci, Reservoir Computing
Approach to Quantum State Measurement, Phys. Rev. X, 11 (2021) 041062.
[19] Pere Mujal, Rodrigo Martínez-Peña, Johannes Nokkala, Jorge García-Beni, Gian Luca
Giorgi, Miguel C. Soriano, Roberta Zambrini, Opportunities in Quantum Reservoir
Computing and Extreme Learning Machines, Advanced Quantum Technologies, 4 (2021)
2100027.
[20] Yudai Suzuki, Qi Gao, Ken C. Pradel, Kenji Yasuoka, Naoki Yamamoto, Natural quantum
reservoir computing for temporal information processing, Scientific Reports, 12 (2022)
1353[21] 有田拳. 卒業論文「超伝導体内の電界の時間変化を用いたリザバーコンピューテ
ィングに関する研究」. 2022 年 3 月.
[22] 松下照男．「磁束ピンニングと電磁現象」第 2 版．2014 年 3 月[23] V.L. Ginzburg and L.D. Landau, Zh. Eksp. Teor. Fiz. 20, 1064 (1950). English translation
in: L. D. Landau, Collected papers (Oxford: Pergamon Press, 1965) p. 546[24] Albert Schmid, A time dependent ginzburg-landau equation and its application to the
problem of resistivity in the mixed state, Physik der kondensierten Materie, 5 (1966) 302–
317.
[25]松下照男．卒論講義ノート- 55 -

[26] Tetsuya Matsuno, Edmund Soji Otabe, and Yasunori Mawatari, Explicit integrators basedon a bipartite lattice and a pair of affine transformations to solve quantum equations withgauge fields, Journal of the Physical Society of Japan, 89 (2020) 054006.
[27] Tetsuya Matsuno, Edmund Soji Otabe, and Yasunori Mawatari, Explicit StructurePreserving Integrators for Dissipative and Conservative Nonlinear Time Dependent
Schrödinger Equations with Gauge Fields, Journal of the Physical Society of Japan, 92(2023) 074004.
[28] T. Matsuno. Link variables for the TDGL equation. September 11, 2015[29] 松野哲也．AFI をはじめよう ver.6．2020 年 1 月[30] 田中剛平，中根了昌, 廣瀬明．リザーバーコンピューティング-時系列パターン認
識のための高速機械学習の理論とハードウェア．森北出版株式会社．2021 年 3 月[31] Amir F Atiya and Alexander G Parlos, New results on recurrent network training: unifying
the algorithms and accelerating convergence, IEEE transactions on neural networks, 11(2000) 697–709.
[32] Herbert Jaeger, Short term memory in echo state networks. gmd-report 152, GMD-German
National Research Institute for Computer Science (2002).

- 56 -

謝辞
本研究に取り組むにあたり、多くの方から多大なご助力を賜りました。
まず、指導教官である小田部荘司教授に御礼申し上げます。研究におきましては進捗が行き詰ったときの相談に乗っていただき、研究に必要になる超伝導現象に関する理論に関しても大変わかりやすくご教授下さいました。また、学会発表など私の研究を発表する機会や、バングラディシュ交流ゼミという英語力を鍛える場を与えていただきました。
次に、有明高専の松野先生に御礼申し上げます。松野先生が書かれた AFI の大変わかりやすい資料のおかげで研究をスムーズに進めることができました。また、私が質問をさせていただいた際にも親切なご回答をいただきました。
次に、ニューロモフィックセンターの宇佐美先生、田中先生に御礼申し上げます。お二方には、にはリザバーコンピューティングに関する部分で大変にお世話になりました。
機械学習などに関して素人である私に対して丁寧に教えてくださいました。また、共同研究までさせていただき、私だけで進めても絶対に得られることのない研究結果を得る
ことができました。
次に、AFI、リザバーコンピューティングの研究テーマに関して先行して研究をされていた上田天馬さんに御礼申し上げます。上田さんが研究資料を細かく残していてくださったおかげで一つひとつしっかりと理解しながら研究を進めることができました。また、私が困ったときには親身に相談に乗ってくださいました。
最後にお世話になった小田部研究室の皆様に深く感謝申し上げます。

- 57 -

研究業績発表
1. 有田拳、上田天馬、小田部荘司、宇佐美雄生、田中啓文、松野 哲也「超電導体内の電界の時間変化を用いたリザバーコンピューティングに関する研究」
2022 年度春季低温工学・超電導学会、タワーホール船堀、令和 4 年 6 月 20 日
2.

K. Arita, T. Ueda, E.S. Otabe, Y. Usami, H. Tanaka, and T. Matsuno. “Possibility of using thepinning phenomenon of superconductors as a reservoir computing”. ISS2022. WINC-AICHI.
December 1, 2022.

3.

有田拳、上田天馬、小田部荘司、宇佐美雄生、田中啓文、松野 哲也「磁束線の運動を利用した超電導リザバーによる音声認識」
2022 年度秋季低温工学・超電導学会、長良川国際会議場、令和 4 年 12 月 9 日

4.

K. Arita, T. Ueda, E.S. Otabe, Y. Usami, H. Tanaka, and T. Matsuno. “Investigation of theeffect of changes in the pinning properties of Type II superconductors on theprediction accuracy of superconducting reservoirs”. The 4th International
Symposium on Neuromorphic AI Hardware. ア ー ト ホ テ ル 小 倉 ニ ュ ー タ ガ ワ ,December 14th, 2022.

- 58 -