Nobisuke

Dekisugi

RAG

autoQAOA

テンソルネットワークツールと基本的な量子ビット演算、特異値分解SVD

Yuichiro Minato

2021/02/13 08:22

#量子ゲート

はじめに

テンソルネットワークという技術を使って量子コンピュータの計算をしてみます。今回は基礎から理解してもらいたいため、簡単なところから始めます。

テンソルネットワーク

テンソルというノードをグラフネットワーク化したもので、ベクトルや行列をグラフのノードとして捉えて計算を行います。

階数

テンソルには階数というものがあり、階数によってテンソルの形が変わります。腕の数を数え、その数によって種類が異なります。順番に見てみます。

階数0、スカラー量

腕がないでノードだけの場合にはスカラー量になります。

無題のプレゼンテーション.jpg

階数1、ベクトル

腕が一本の場合にはベクトルになります。

無題のプレゼンテーション (1).jpg

v = \begin{bmatrix} v_0\\ v_1\\ \vdots\\ v_i \end{bmatrix}

階数2、行列

腕が二本の場合には行列になります。腕の両側の次元が行列の縦横の次元に対応します。

無題のプレゼンテーション (2).jpg

M = \begin{bmatrix} M_{00} & M_{10} & \cdots & M_{j0}\\ M_{01} & & &\\ \vdots\\ M_{0i} & & & M_{ji} \end{bmatrix}

階数3,4、それぞれの階数のテンソル

腕が三本、四本と増えるにつれてテンソルの次元が増えていきます。

無題のプレゼンテーション (3).jpg

テンソルの縮約

それぞれのテンソルを繋ぎますが、繋ぐ腕の次元が等しい場合、縮約と言ってまとめることができます。

無題のプレゼンテーション (4).jpg

行列同士は行列に、行列とベクトルはベクトルになります。

無題のプレゼンテーション.jpg

その他次元の高いテンソルも同様に腕の数を見ます。

テンソルの分解（特異値分解）

テンソルは分解もできます。特異値分解(Singular Value Decomposition)SVDを使います。対角行列１つと直交行列２つに分解できます。

無題のプレゼンテーション (1).jpg

次はツールの使い方を見てみます。

google tensornetwork

こちらは勉強会でもおなじみのgoogleのtensornetworkです。

https://github.com/google/TensorNetwork

今回はさらに計算を速度向上させるようにバックエンドについてちょっと見てみましょう。

インストールは、

pip3 install tensornetwork

基本的な使い方は、

import numpy as np
import tensornetwork as tn

a = tn.Node(np.ones((10,))) 
b = tn.Node(np.ones((10,)))
edge = a[0] ^ b[0]
final_node = tn.contract(edge)
print(final_node.tensor)

こんな感じでノードを作って縮約を取ります。今回さらに強化してみます。

!pip install tensornetwork jax jaxlib

import numpy as np
import jax
import tensornetwork as tn

今回は二つのバックエンドの速度を比較してみます。

def calculate_abc_trace(a, b, c):
  an = tn.Node(a)
  bn = tn.Node(b)
  cn = tn.Node(c)
  an[1] ^ bn[0]
  bn[1] ^ cn[0]
  cn[1] ^ an[0]
  return (an @ bn @ cn).tensor

a = np.ones((4096, 4096))
b = np.ones((4096, 4096))
c = np.ones((4096, 4096))

tn.set_default_backend("numpy")
print("Numpy Backend")
%timeit calculate_abc_trace(a, b, c)
tn.set_default_backend("jax")

print("JAX Backend")
%timeit np.array(calculate_abc_trace(a, b, c))

GPUモードにするの忘れてました。

Numpy Backend
1 loop, best of 3: 4.5 s per loop
JAX Backend
/usr/local/lib/python3.6/dist-packages/jax/lib/xla_bridge.py:119: UserWarning: No GPU/TPU found, falling back to CPU.
  warnings.warn('No GPU/TPU found, falling back to CPU.')
1 loop, best of 3: 2.93 s per loop

それでもちょっとだけ早いです。GPUでやってみます。

Numpy Backend
1 loop, best of 3: 3.6 s per loop
JAX Backend
The slowest run took 28.18 times longer than the fastest. This could mean that an intermediate result is being cached.
1 loop, best of 3: 104 ms per loop

かなり早くなりました。

quimb

こちらはグラフを見てみます。

pip install quimb autoray

%matplotlib inline
import quimb as qu
import quimb.tensor as qtn

# 10 qubits and tag the initial wavefunction tensors
circ = qtn.Circuit(N=10, tags='PSI0')

# initial layer of hadamards
for i in range(10):
    circ.apply_gate('H', i, gate_round=0)

# 8 rounds of entangling gates
for r in range(1, 9):
    # even pairs
    for i in range(0, 10, 2):
        circ.apply_gate('CNOT', i, i + 1, gate_round=r)

    # Y-rotations
    for i in range(10):
        circ.apply_gate('RY', 1.234, i, gate_round=r)

    # odd pairs
    for i in range(1, 9, 2):
        circ.apply_gate('CZ', i, i + 1, gate_round=r)

# final layer of hadamards
for i in range(10):
    circ.apply_gate('H', i, gate_round=r + 1)

circ
circ.psi.graph(color=['PSI0', 'H', 'CNOT', 'RY', 'CZ'])

ちょっと複雑なゲートを実装してみました。

ダウンロード.png

色によってゲートが規定されていますが、いい感じですね。

quflex

https://github.com/ngnrsaa/qflex

こちらはNASAですね。いうほどテンソルネットワーク使ってないみたいですが、一応紹介しておきました。 ##Google Tensornetworkを使ってみる Googleからツールが出ています。

pip install tensornetwork

基本操作を確認してみましょう。

import numpy as np
import tensornetwork as tn

#ベクトルを作る
a = tn.Node(np.ones((10,))) 
b = tn.Node(np.ones((10,)))

#ベクトル同士を繋ぐ
edge = a[0] ^ b[0] # Equal to tn.connect(a[0], b[0])

#縮約
final_node = tn.contract(edge)

#結果の表示
print(final_node.tensor) # Should print 10.0

このようにテンソルを縮約してベクトル２つからスカラー量がもとまりました。

量子ゲート計算をしてみる

テンソルネットワークは量子コンピュータの計算のためにあるわけではありませんが、量子コンピュータではベクトルとしての「状態ベクトル」とユニタリ行列としての「量子ゲート」操作があるため、全ての操作をテンソルネットワークで表現できます。

量子ビット

量子ビットは|0>は通常状態ベクトルの[1,0]の縦ベクトルで表現されます。

import numpy as np
import tensornetwork as tn

q0 = tn.Node(np.array([1,0]))

早速これに量子ゲートをかけてみます。量子ゲートはXゲートやHゲートを作用させてみます。

X = tn.Node(np.array([[0,1],[1,0]]))
H = tn.Node(np.array([[1,1],[1,-1]])/np.sqrt(2))

q0 = tn.Node(np.array([1,0])) 

edge1 = q0[0] ^ X[0]
print(tn.contract(edge1).tensor)

edge2 = q0[0] ^ H[0]
print(tn.contract(edge2).tensor)

答えは、それぞれ

[0 1]
[0.70710678 0.70710678]

と、それぞれXゲート、Hゲートの適用になりました。今度はアダマールからZゲート行きます。

q = tn.Node(np.array([1,0])) 
Z = tn.Node(np.array([[1,0],[0,-1]]))

q[0] ^ H[0]
H[1] ^ Z[0]

a = q@H
b = a@Z

print(b.tensor)

次に再度HからXゲート行ってみます。

H = tn.Node(np.array([[1,1],[1,-1]])/np.sqrt(2))
X = tn.Node(np.array([[0,1],[1,0]]))

b[0]^H[0]
H[1]^X[0]

c = b@H
d = c@X

print(c.tensor)
print(d.tensor)

それぞれ、

[0. 1.]
[1. 0.]

一周回ってきました。

CNOT

CNOT回路を実装しますが、

q0 *--  --*--
          *   CX
q1 *--  --*--

q0とq1のベクトルを準備して、CXを階数４のテンソルとして準備します。つなげて、

q0 *----*--
        * 
q1 *----*--

順番に縮約を取ります。

!pip install tensornetwork

import numpy as np
import tensornetwork as tn

q0 = tn.Node(np.array([1,0])) 
q1 = tn.Node(np.array([1,0]))

X = tn.Node(np.array([[0,1],[1,0]]))
CX = tn.Node(np.array([[1,0,0,0],[0,1,0,0],[0,0,0,1],[0,0,1,0]]).reshape((2,2,2,2)))

q0[0] ^ CX[0]
q1[0] ^ CX[1]

cx1 = q0@CX
cx2 = q1@cx1

print(cx2.tensor.reshape(4))

#=>[1 0 0 0]

こちらは何も起きず、|00>のままです。

コントロールビットを1に

次にコントロールビットを１にしてみます。

import numpy as np
import tensornetwork as tn

q0 = tn.Node(np.array([1,0])) 
q1 = tn.Node(np.array([1,0]))

X = tn.Node(np.array([[0,1],[1,0]]))
CX = tn.Node(np.array([[1,0,0,0],[0,1,0,0],[0,0,0,1],[0,0,1,0]]).reshape((2,2,2,2)))

q0[0] ^ X[0]
X[1] ^ CX[0]
q1[0] ^ CX[1]

x1 = q0@X
cx1 = x1@CX
cx2 = q1@cx1

print(cx2.tensor.reshape(4))

#=>[0 0 0 1]

波動関数が下記かわり|11>となりました。

最後にもつれ

量子もつれを作ってみます。

*--H--*--
      |
*-----X--

q0 = tn.Node(np.array([1,0])) 
q1 = tn.Node(np.array([1,0]))

H = tn.Node(np.array([[1,1],[1,-1]])/np.sqrt(2))
CX = tn.Node(np.array([[1,0,0,0],[0,1,0,0],[0,0,0,1],[0,0,1,0]]).reshape((2,2,2,2)))

q0[0] ^ H[0]
H[1] ^ CX[0]
q1[0] ^ CX[1]

h1 = q0@H
cx1 = h1@CX
cx2 = q1@cx1

print(cx2.tensor.reshape(4))

#=>[0.70710678 0.         0.         0.70710678]

MPS（行列積状態）

状態ベクトルの表現を行列積と呼ばれる状態に持っていくことで計算がしやすくなります。通常量子ビットを使って量子状態の表現は $2^N$ の計算量が必要ですが、階数３のテンソルを使って量子状態を表すことで、 $2Nm^2$ 程度で計算が済みます。 $m$ はテンソル間の結合次元が対応します。

無題のプレゼンテーション.jpg

##結合次元の圧縮 SVDした時に対角行列sをみることで、次元を圧縮できる。

U
\begin{bmatrix}
a&&&\\
&a&&\\
&&0&\\
&&&0
\end{bmatrix}
V

SVDのコード

こちらです。バックエンドによってちょっとずつ使い方が違います。

import tensornetwork as tn
node = tn.Node(np.array([[1.,2.],[3.,4.]]), backend="jax")
U, s, V, trun_error = tn.split_node_full_svd(node, [node[0]], [node[1]])

U.tensor
DeviceArray([[-0.40455365, -0.91451436],
             [-0.9145144 ,  0.4045536 ]], dtype=float32)

import tensornetwork as tn
node = tn.Node(np.array([[1.,2.],[3.,4.]]),backend="numpy")
U, s, V, trun_error = tn.split_node_full_svd(node, [node[0]], [node[1]])
s.tensor

array([[5.4649857 , 0.        ],
       [0.        , 0.36596619]])

import tensornetwork as tn
import tensorflow as tf

node = tn.Node(np.array([[1.0,2.0],[3.0,4.0]]), backend="tensorflow")
U, s, V, trun_error = tn.split_node_full_svd(node, [node[0]], [node[1]])

sess= tf.Session()

print(V.tensor.eval(session=sess))
[[ 0.57604844  0.81741556]
 [ 0.81741556 -0.57604844]]

print(s.tensor.eval(session=sess))
[[5.4649857  0.        ]
 [0.         0.36596619]]

print( (U @ s @ V).tensor.eval(session=sess))
[[1. 2.]
 [3. 4.]]

split node

こちらは

import tensornetwork as tn
import numpy as np
node = tn.Node(np.array([[1.,2.],[3.,4.]]),backend="numpy")
U, V, trun_error = tn.split_node(node, [node[0]], [node[1]])

U.tensor

sの平方根がとられてUとVに自動導入されます。使いどころあったりなかったり。

numpyでもやってみます

numpyでも分解してみます。最初はnumpyを使ってみます。

import numpy as np
A = np.array([[1,2],[3,4]])
U, s, V = np.linalg.svd(A, full_matrices=True)

これで、

U
#=>array([[-0.40455358, -0.9145143 ],
       [-0.9145143 ,  0.40455358]])

s
#=>array([5.4649857 , 0.36596619])

V
#array([[-0.57604844, -0.81741556],
       [ 0.81741556, -0.57604844]])

注意点はsが一次元の配列で出てくるので、行列として表現するには、

np.diag(s)

array([[5.4649857 , 0.        ],
       [0.        , 0.36596619]])

として確認できます。確かめは、元に戻せるかでできます。

U@np.diag(s)@V

#=>array([[1., 2.],
       [3., 4.]])

おまけ(numpyの多次元配列)

テンソルを扱おうとすると行列の次元を扱う必要があります。ここではちょっとみてみます。numpyを使います。

import numpy as np

１次元配列

print(np.arange(3).reshape((3)))


[0 1 2]

ベクトルになりました。

2次元配列

print(np.arange(9).reshape((3,3)))


[[0 1 2]
 [3 4 5]
 [6 7 8]]

行列です

3次元配列

print(np.arange(27).reshape((3,3,3)))

[[[ 0  1  2]
  [ 3  4  5]
  [ 6  7  8]]

 [[ 9 10 11]
  [12 13 14]
  [15 16 17]]

 [[18 19 20]
  [21 22 23]
  [24 25 26]]]

2次元の配列が並んでいます。イメージとしては奥行きが出た感じでしょうか。

４次元配列

print(np.arange(81).reshape((3,3,3,3)))

[[[[ 0  1  2]
   [ 3  4  5]
   [ 6  7  8]]

  [[ 9 10 11]
   [12 13 14]
   [15 16 17]]

  [[18 19 20]
   [21 22 23]
   [24 25 26]]]


 [[[27 28 29]
   [30 31 32]
   [33 34 35]]

  [[36 37 38]
   [39 40 41]
   [42 43 44]]

  [[45 46 47]
   [48 49 50]
   [51 52 53]]]


 [[[54 55 56]
   [57 58 59]
   [60 61 62]]

  [[63 64 65]
   [66 67 68]
   [69 70 71]]

  [[72 73 74]
   [75 76 77]
   [78 79 80]]]]

上記の三次元の奥行きを持った行列の塊がさらにできました。テンソルネットワークでは当面４くらいまで覚えればなんかできそうな気がします。

テンソルの設定

googleのtensornetworkを使って設定します。

import numpy as np
import tensornetwork as tn

q0 = tn.Node(np.arange(81).reshape((3,3,3,3))) 
print(q0.tensor)

無事設定できました。

[[[[ 0  1  2]
   [ 3  4  5]
   [ 6  7  8]]

  [[ 9 10 11]
   [12 13 14]
   [15 16 17]]

  [[18 19 20]
   [21 22 23]
   [24 25 26]]]


 [[[27 28 29]
   [30 31 32]
   [33 34 35]]

  [[36 37 38]
   [39 40 41]
   [42 43 44]]

  [[45 46 47]
   [48 49 50]
   [51 52 53]]]


 [[[54 55 56]
   [57 58 59]
   [60 61 62]]

  [[63 64 65]
   [66 67 68]
   [69 70 71]]

  [[72 73 74]
   [75 76 77]
   [78 79 80]]]]

これで色々捗りそうです。

テンソルネットワークツールと基本的な量子ビット演算、特異値分解SVD

Yuichiro Minato

はじめに.css-1yssqf7{opacity:0.6;border:0;border-color:inherit;border-style:solid;border-bottom-width:1px;width:100%;margin-top:var(--chakra-space-2);}

テンソルネットワーク

階数

階数0、スカラー量

階数1、ベクトル

階数2、行列

階数3,4、それぞれの階数のテンソル

テンソルの縮約

テンソルの分解（特異値分解）

google tensornetwork

quimb

quflex

量子ゲート計算をしてみる

量子ビット

CNOT

コントロールビットを1に

最後にもつれ

MPS（行列積状態）

SVDのコード

split node

numpyでもやってみます

おまけ(numpyの多次元配列)

１次元配列

2次元配列

3次元配列

４次元配列

テンソルの設定

はじめに