量子状態トモグラフィのチュートリアル(Qiskit)をやってみる

この記事は量子コンピュータ Advent Calendar 2020の12/9ぶんです。
qiita.com


こんにちは。Kumaです。
この記事では、量子状態を測定で特定する方法を考えたいと思います。
使う手法は「量子状態トモグラフィ」と呼ばれるものです。

量子状態トモグラフィ

量子状態トモグラフィ　を説明するには、量子状態の密度演算子とそのパウリ行列展開を考える必要があります。
electrodynamics.hatenablog.com
qiita.com


量子状態は一般に複素ベクトル | \psi \rangle で表現できますが、密度演算子表示というものもあります。
密度演算子は、 \rho = | \psi \rangle \langle \psi |　で定義されます。（ベクトルの外積なので、行列になります）
密度演算子はエルミート演算子であって、パウリ行列 \sigma_{x}, \sigma_{y}, \sigma_{z}を基底として展開できます。
その座標成分（展開したときの係数）は必ず実数となります。
注目したいのは、これらの係数が観測可能な物理量となっていることです。
たとえば1量子ビットの場合、係数は（定数を除けば）3つあるのですが、それぞれ x,y,z方向のスピンの期待値 \langle \psi | \sigma_{x,y,z} | \psi \rangleに対応します。
スピンの期待値を測ることで、係数が求まり、密度演算子が判明します。

ここで重要な注意ですが、 x,y,z方向のスピンの期待値は、同時には測定できません。
なぜならば、 x,y,z方向のスピンを表す演算子は可換ではないからです。
また、期待値を求めなければならないので、状態の生成と測定を繰り返す必要があります。（測定のつど、状態は射影されて変わってしまうので、再度生成が必要となる）
一般には、雑音もある実回路で毎回同じ状態を用意できるとは限りませんので、これも注意が必要なポイントです。

Qiskitによるチュートリアル

Qiskitには量子状態トモグラフィのチュートリアルもあります。
qiskit.org

これを実行してみます。まず、pip install qiskit , pip install pylatexenc の後に、importしていきます。

f:id:phymath1991:20201208031434p:plain

次に、推定してみたい量子状態を定義します。
f:id:phymath1991:20201208031556p:plain
f:id:phymath1991:20201208031623p:plain

量子状態トモグラフィをかける回路を作ります。
f:id:phymath1991:20201208031713p:plain
ここでわざわざ冗長な6qubitの回路を用いる意味がよくわかっていません。。

基底（パウリ行列）の期待値を測定します。
そして、StateTomographyFitterに引き渡します。
f:id:phymath1991:20201208031800p:plain

密度演算子の推定を実施します。
qiskitのhelpをみると、method='lstsq'は、「最小二乗規範で推定します」ということのようです。
このあたりもよくわかっていませんが。
推定した状態と、答えの状態のフィデリティを計算します。フィデリティは1に近いほどよいです。
f:id:phymath1991:20201208032017p:plain

だいたい1に近いので、推定できていると言えますね。
ちなみに、期待値測定の結果は、実行の度変わります。（無限回測定しているわけではないので、ゆらぎがある）
よって、フィデリティも微妙に変わります。
f:id:phymath1991:20201208032335p:plain
f:id:phymath1991:20201208032354p:plain

実機でもやってみたいですね。

私の量子コンピュータ史

この記事は量子コンピュータ Advent Calendar 2020の12/1ぶんです。
qiita.com


こんにちは。Kumaです。
この記事では「私の量子コンピュータ史」を書いてみたいと思います。*1

量子コンピュータに出会うきっかけ

私が量子コンピュータの勉強やプログラミングをやろうと思ったきっかけは、
趣味で参加した　日曜数学会　というイベントでした。
このイベントの第12回目が2018年6月に開かれました。
そこで、私は衝撃を受けるのです！
www.nicovideo.jp

こちらの発表は @snuffkin さんによるもので、量子フーリエ変換の解説になります。
まず、量子コンピュータが個人でプログラミングできるということに大変驚きました。
痺れました。
発表後の懇親会タイムで、すぐに@snuffkin　さんに声をかけてお話を伺いました。

@snuffkin　さんによると、もし興味があるのであれば、紹介したい人がいる　とのことでした。
そこで出会ったのが　@kyamaz　さんでした。
@kyamaz　さんは、OpenQL　という量子コンピュータのライブラリ等のオープン開発プロジェクトを主催している方です。
openql.connpass.com

OpenQL　は勉強会もやっており、ぜひ一度見に来てくださいとお誘いを頂きました。

この出会いは非常に幸運だったと思います。
量子コンピュータのプログラミングをやりたいと思ったその場で、勉強するきっかけを作ることが出来ました。
大学で量子力学を独学で身につけた*2自分にとって、これはその腕を振るうまたとないチャンスだと思いました。

OpenQL 勉強会

OpenQL 勉強会　では、@kyamaz　さんが基本的な部分をイントロしてくださるので、ここで概要が一気に掴めました。
特に、量子コンピュータの理想と現状の差は大きいこと、多様な分野の専門家が参加している（参加する必要がある）こともわかりました。
OpenQL勉強会では、聴講だけでなく、ライトニングトーク（LT)も何度かやらせていただきました。
勉強するモチベーションとして、「まとめてアウトプットする」は必須な作業なので、OpenQLの存在は非常に大きかったです。
この勉強会では、ゲートの他にもアニーリングなども聴講できるので、とにかく勉強になりました。
アニーリングをやっている人に色々聞けたことは大きかったですね。
また、　@Ayumu_walker さんがQiskit adovocateや Qiskit camp の話をしていて、これも大変おもしろかったです。*3
qiskit.org

他にも、コミケに参加して量子コンピュータの解説本を買いあさりました。
なお、@snuffkin　さんは本も書かれています。なんでも出来るなこの人！
booth.pm

日曜数学会へのフィードバック

量子コンピュータを少しずつ勉強していくと、基本的な部分は理解できるようになりました。
しかし、一般人が「量子コンピュータをやっている人に量子コンピュータの話を聞く機会」は、どれぐらいあるでしょうか？
ほとんどないのでは？と思った私は、「自分が理解した範囲の量子コンピュータの話を、原点である　数学会　にフィードバックしたらみんなどんな反応をするのだろう？」と考えました。
まさに　@snuffkin さん　をリスペクトすることにしたのです。

実はわたしは、友人数人と　日曜数学会（東京）の分会である　関西日曜数学友の会　を主催しております。
ちょうど第3回の招待講演者を探していたところだったので、
私の一存（ゴリ押し）で、OpenQLの　@kyamazさんを招待講演者として召喚しました！
職権乱用！

techplay.jp

私も、量子コンピュータの基礎的な部分をLTしました。
Q&Aも大いに盛り上がり、量子コンピュータへの関心の高さ（しかし敷居が高いのでなかなか情報が得られないもどかしさ）が伺えました。
感触を掴んだ私は、本腰を入れて量子コンピュータの勉強をやろうと考え、ニールセンを購入しました。*4
www.amazon.co.jp


最近の日曜数学会（東京）でも、量子ゲートの話をしました。*5
非常に基本的な部分ではありますが。
www.nicovideo.jp
electrodynamics.hatenablog.com

今後も、工夫をこらして量子コンピュータの話題をねじ込んでいこうと思います。

今の話

今も基礎的な部分の勉強を続けています。
最近はqiskitなどを使いながら、実機実行も含めてちまちまやっています。
本業の傍らなので進みは遅いですが、気長にやっていくつもりです。
なにか量子コンピュータの進展に貢献できる・・・かはわかりませんが、まずは楽しむことを頑張ろうと思います。

以上、私の量子コンピュータ史でした。
良いお年を。

*1:量子コンピュータはいろいろなバックグラウンドを持つ（かつ、新しいものに興味がある）方が参加されていると思うので、 他の方の話も聞いてみたいですね。

*2:身につけた、といっても、基本が分かる程度ではありますが。

*3:"時を巻き戻すアルゴリズム"" というキャッチーなタイトルのLTもありました。ネーミングすこ。

*4:このあたりで、量子コンピュータでは量子情報理論、特に誤り訂正符号が非常に重要なことを知りました。たまたま本業で情報理論や誤り訂正符号を少しだけやったことのあった私は、運命のめぐり合わせを感じたのです。勉強することは山程ありますが。

*5:発表中に、娘が泣いていたので音声が入ってしまっています。子育てパパでも量子プログラミングは出来るよ！

多量子ビット状態の可視化(多次元ブロッホ球)

この記事は量子コンピュータ Advent Calendar 2020の12/3ぶんです。
qiita.com


こんにちは。Kumaです。
この記事では量子ビット状態の可視化を考えてみます。
可視化、というと、1量子ビット状態ではブロッホ球なるものが出てきます。
しかし、2量子ビット状態になると、ブロッホ球の話は全然出てこなくなります。
なぜでしょうか？？

実は、1量子ビットの可視化は容易なのですが、2量子ビットに増えるだけで一気に難しくなります。

道具の準備

天下り的ですが、必要な道具を先に導入します。

2x2エルミート行列 Hは、パウリ行列 \sigma_{1}, \sigma_{2}, \sigma_{3} と呼ばれる特別なエルミート行列と、単位行列 I の線形結合で書くことができます。
H = \frac{1}{2} ( h_{0}I + h_{1}\sigma_{1} + h_{2}\sigma_{2} + h_{3}\sigma_{3} )
ここで h_{0}, h_{1}, h_{2}, h_{3} はすべて実数であって次の式を満たします。
h_{0} = \rm{Trace} (IH), h_{1} = \rm{Trace}(\sigma_{1}H), h_{2} = \rm{Trace}(\sigma_{2}H), h_{3} = trace(\sigma_{3}H)
ちなみに、Traceというのは行列の対角成分全ての和を取る演算です。

f:id:phymath1991:20181209203724p:plain:w500
エルミート演算子のパウリ行列展開
これは以前の記事で証明しています。
electrodynamics.hatenablog.com


他方、二成分の複素ベクトル | \psi \rangle

| \psi \rangle = \begin{pmatrix} \alpha \\ \beta \end{pmatrix}

に対して

| \psi \rangle \langle \psi | = \begin{pmatrix} \alpha \alpha^{\star} & \alpha \beta^{\star} \\ \beta \alpha^{\star} & \beta \beta^{\star} \end{pmatrix}
を | \psi \rangleの密度演算子といい、記号 \rho で書きます。（このブログでは　射影演算子　という呼び方もしていますが、同じです）
密度演算子は明らかにエルミート行列なので、パウリ行列展開が可能です。

1量子ビット状態は1つのブロッホ球

密度演算子のパウリ行列展開における係数(”成分”)である実数 \psi_{0,1,2,3}を考えます。
ここで のノルム |\alpha|^{2}+|\beta|^{2}=1であるとします。
すると | \psi \rangleによらず \psi_{0} = 1であることがわかるので、 \psi_{0}は定数です。自由度ではないので無視しましょう。
結局、密度演算子を考えることでノルム１の２成分複素ベクトル | \psi \rangleは 実数の三次元ベクトル
\vec{\psi} = \begin{pmatrix} \psi_{1} \\ \psi_{2} \\ \psi_{3} \\ \end{pmatrix}
に対応しました。
さらにこれらは独立ではありません。
それを確認するために、成分表示を定義から求めてみると
\psi_{1} = |\alpha|^{2} - |\beta|^{2} \\ \psi_{2} = \alpha^{\star}\beta + \beta^{ \star}\alpha = 2Re[\alpha \beta^{\star} ] \\ \psi_{3} = -j(\alpha^{\star}\beta -\beta^{ \star}\alpha) = -2Im[\alpha \beta^{\star} ] \\
となります。上記からも
\psi_{1}^{2}+\psi_{2}^{2}+\psi_{3}^{2} = 1 という関係式を満たすことがわかります。
すなわち、 \vec{\psi} は三次元空間における半径１の球面と１：１に対応します。
*1

この球をブロッホ球、 \vec{\psi} をブロッホベクトルといいます。*2

f:id:phymath1991:20181209204622p:plain:w500
射影演算子の”座標成分”の球面表示

2量子ビット状態は2つのブロッホ球？

さて、2量子ビットに拡張してみましょう。1量子ビット状態が２つあって、それぞれ | \psi \rangleと | \phi \rangleだったとします。
これを2量子ビット状態だとみなすとき、その量子状態は、 \otimes という記号を用いて、
| \psi \rangle \otimes | \psi \rangleと書けます。
よくわかりませんが、成分で定義すると
| \psi \rangle \otimes | \psi \rangle = \begin{pmatrix} \alpha \gamma \\ \alpha \delta \\ \beta \gamma \\ \beta \delta \end{pmatrix}
ということです。ここで \alpha, \beta は | \psi \rangleの成分で、 \gamma, \delta は | \phi \rangleの成分です。
2量子ビット状態は、ノルム1の4成分複素数ベクトルで書るということですね。
単純に考えれば、成分の数が2倍になっているので、ブロッホ球を2個並べれば表現できそうですね。実際これは正しいです。お疲れ様でした。
と、いうわけにはいきません。

2量子ビット状態をちゃんと考えてみる

なぜブロッホ球を2個並べてもだめなのでしょうか？
ブロッホ球を2個並べて表示できるのは、1量子ビット目と2量子ビット目が無相関である時だけだからです。
実際、我々は意図せずに無相関であることを仮定してしまっています。

我々の先の数式展開に紐づけてみると、以下が問題なわけです。

1量子ビット状態が２つあって、それぞれ | \psi \rangleと | \phi \rangleだったとします。
これを2量子ビット状態だとみなすとき、その量子状態は、 \otimes という記号を用いて、
| \psi \rangle \otimes | \psi \rangleと書けます。

ここがまずいのです。
このような状態は積状態と呼ばれます。積状態は、ただ1量子ビット状態を２つもってきただけなので、互いに相関はありません。
しかし積状態は量子状態全てを尽くしていません。
例えば、任意の積状態の和もまた新たな量子状態として認めなければなりません。
例えば、2量子ビット積状態を2つ考えて、その重ね合わせ | \psi \rangle \otimes | \phi \rangle + | \lambda \rangle \otimes | \eta \rangle も2量子ビットの状態として存在しえます。
しかしこれは一般に　 | \dot \rangle \otimes | \dot \rangle　のような積状態に分解が出来ません。
特に極端なものとして、 | \psi \rangle = \frac{1}{\sqrt{2}} ( | 00 \rangle + | 11 \rangle ) という状態（エンタングルメント）があります。
これは第一量子ビットと第二量子ビットが常に同じ値を取るように束縛がかかった状態です。（もちろん、積状態で書けません）
このように、量子状態では、互いの状態が相関することがありえます。例えば1量子ビット目がなにか変化したときに、2量子ビット目が巻き込まれて変わってしまうことがあるのです。

よって、「積状態であれば、(独立な）ブロッホ球を2個並べて可視化できる」のですが、一般の2量子ビット状態では、そのようなことはできません。
”相関する２つのブロッホ球”　を考えなければなりません。

2量子ビット状態のブロッホベクトル

では、2量子ビット状態のブロッホベクトルはどうあるべきなのでしょうか？
以下の論文を参考にして考えました。
[1602.01548] Entangled Bloch Spheres: Bloch Matrix and Two Qubit State Space

まず、2量子ビット状態と、その密度演算子 \rhoを考えます。
性懲りもなく、密度演算子のパウリ行列展開を試みます。
実は基底として \sigma_{\mu} \otimes \sigma_{\nu}、 \mu , \nu = 0,1,2,3 を取ると、展開できます。
すなわち、 \rho = \sum_{\mu , \nu = 0,1,2,3} \rho_{\mu \nu} \sigma_{\mu} \otimes \sigma_{\nu} ここで \rho_{\mu \nu}は全て実数です。
すると、ブロッホベクトルならぬ、4x4の”ブロッホ行列” \rho_{\mu \nu}が状態を決定づけそうです。
ブロッホ行列を少し調べてみます。成分計算をしていくと、ブロッホ行列の中身は次のように分類できます。
\rho = \begin{pmatrix} 1 && \vec{v}^{\dagger} \\ \vec{u} && R \end{pmatrix}
ここで \vec{u} は第一量子ビットに対する（ローカル）ブロッホベクトル、 \vec{v}は第二量子ビットに対する（ローカル）ブロッホベクトルです。
どちらも3次元の実数ベクトルです。
Rは3x3の行列で、相関行列と呼ばれます。

積状態のブロッホベクトル

一番簡単な状態として、積状態を考えます。このとき、Rは \vec{u} , \vec{v}だけで完全に決まります。具体的には　 R = \vec{u} \vec{v}^{\dagger} となります。（ベクトル外積であることに注意）
よって、ブロッホ行列は単に各量子ビットに対する（独立な）ブロッホベクトルで完全に指定されますので、独立な2つのブロッホ球で表現できます。

しかし一般の状態では、Rを考慮する必要があることは明らかです。よって、独立な2つのブロッホ球では記述できないわけです。

一般の状態のブロッホベクトル

ここで、4x4のブロッホ行列を4x4=16個の実数成分を持つ1つのベクトルとみてみます。
更に、ブロッホ行列の(0,0)成分は常に定数1だったので、無視します。（15次元のベクトルになります）
このやり方は、以下の（おそらく多次元ブロッホ球で最も有名な）論文にあります。*3
[quant-ph/0301152] The Bloch Vector for N-Level Systems

このとき、ブロッホベクトルはノルムが常に1でなければならないことがわかります。
1量子ビット系のときと同じですね。
すなわち、ブロッホベクトルは15次元球面に存在します。*4

では、15次元球面の任意の点は、意味ある量子状態 | \phi \rangle でしょうか？
1量子ビット状態のときはそうだったので、2量子ビットでも当然そうであるように思います。
しかし、これは間違いです！

詳細は先の論文を見て頂くとよいのですが、15次元球面の一部の領域だけが意味ある量子状態 | \phi \rangle になっています。*5
これは不思議なことでは有りません。そもそも2量子ビットの状態は、複素数4つからなる複素ベクトル | \phi \rangle だったので、
実数にして4x2 = 8 次元しかありません。*6
そのため、15次元球面の全てを埋め尽くせそうにないことがわかると思います。
ブロッホ球は穴だらけなのです。

例えば、二次元射影してやると、以下のようになります。*7
もし球面全体であれば、その射影としては円が浮かび上がるはずですが、そうなりません。

f:id:phymath1991:20201129233313p:plain
15次元ブロッホ球（内部も含んだ）の二次元射影。円にならない。

二量子ビットのブロッホベクトルでも、全然自明な幾何学を持ちません。
難しいですね！
なぜ、ブロッホベクトルは1量子ビット系の時しか登場しないのか、おわかりいただけたでしょうか・・・。

まとめ

・1量子ビット状態は、(三次元）ブロッホ球の表面の点として表現できる。
・2量子ビット状態になると、15次元球面の点として表現できる。（次元が一気にあがる！）
しかも、その球面の一部の点だけが意味ある量子状態を持つ。（ブロッホ球は穴だらけ）
・積状態に限れば、量子ビット間は無相関なので(三次元）ブロッホ球を複数並べればいい。
-------------------------------------------------------------------

*1:ちなみに、 \vec{\psi} の成分表示を一度忘れて、球の内部に相当する量子状態 \vec{\psi} も考えることはできます。この場合、もともとの複素２成分ベクトル | \psi \rangleはどんなものになっているでしょうか？　ふつうに考えると、対応する | \psi \rangleは存在しないのですが、混合状態というものを認めると矛盾なく説明ができるようになります。ここでは混合状態に立ち入ることはしません。

*2:結晶格子の分野でもブロッホというワードが出てきますが、無関係です

*3:ちなみに、この論文は途中計算に大きな間違いがありますので注意してください。

*4:なお、一般には N^{2}-1次元球面になります。いまは N=4です。

*5:最も極端な例として、 \vec{S}がある量子状態に対応するブロッホベクトルであるとき、その負ベクトル -\vec{S}に対応する量子状態は存在しません。これは、ブロッホベクトル同士の内積と量子状態の関係を使うと示せます。

*6:ノルムが1であること、グローバル位相が無視されることなどを要請すると、もっと減ります

*7:ランダムに量子状態を発生させて、そのブロッホベクトルを計算してプロットしています

KumaのLTまとめ "今さら聞けない 量子ゲートコンピュータ"　(第19回　日曜数学会 )

こんにちは。Kumaです。
第19回の日曜数学会で発表した
「今さら聞けない 量子ゲートコンピュータ」
について、内容の振り返りとQ/Aで返しきれなかったものをまとめます。

Title

f:id:phymath1991:20201031160122p:plain

モチベーション

f:id:phymath1991:20201031160147p:plain

Agenda

f:id:phymath1991:20201031160217p:plain

古典状態と量子状態

f:id:phymath1991:20201031160239p:plain

補足*1

量子計算のイメージ

f:id:phymath1991:20201031160326p:plain
f:id:phymath1991:20201031160332p:plain

Groverのアルゴリズム

f:id:phymath1991:20201031160342p:plain
f:id:phymath1991:20201031160434p:plain
f:id:phymath1991:20201031160600p:plain

補足*2
補足2*3

実機でやる

f:id:phymath1991:20201031160631p:plain
f:id:phymath1991:20201031160642p:plain
f:id:phymath1991:20201031160650p:plain

まとめ

f:id:phymath1991:20201031160658p:plain

Q&A

コメントへの回答です

  • 測定は良いランダムなんですか？

→　量子現象なので、完全にランダムです。（この世で最もランダムです。疑似乱数ではなく、真の乱数です。）
物理乱数のランダム性は、熱雑音などにも見られますので、量子測定特有のものというわけではありません。
ただし、量子コンピュータの実機が真にランダムであるかは（意図しない偏りを持つ雑音などが混入する余地があります）わかりません。
量子コンピュータを真の乱数生成器としてとらえたアプリケーションもあります。例えば”量子星占い”があります。
https://www.quantumcomputer.tokyo/horoscope.html

→　Yes。無料で、IBMQにTwitterアカウント等でログインすればGUIが使えます。
https://quantum-computing.ibm.com/
GUI上では、ドラッグアンドドロップでゲート回路の作成や実機の実行が可能です。
しかし5qubitまでしか使えない等の制約があります。
本格的に遊びたい場合は、ぜひプログラミング言語経由での開発をおすすめします。
なお、Google colab. サービスを使うことで、GoogleアカウントさえあればWebブラウザ上でプログラミングの実行が出来ます。

  • あんまり知らないけど組合せ論的なの早い感じなの

→　組み合わせ最適化は良い応用例と考えられています。例えばQAAやQAOAなどのアルゴリズムがあります。
特にQAOAは雑音がある量子コンピュータ(NISQ)でも動くように工夫されたものです。
ただ、問題の複雑度やQAOAのゲートの深さを上げていくと、正解率がどんどん下がってしまうことが知られています。
雑音のある量子コンピュータはスケーラビリティに欠ける印象があります。
しかしある程度のスケーラビリティがないと、古典スパコンで解ける範囲を出られないという問題があります。

→　雑音ありでも動くかも？という観点で、非常に注目されてます。ある程度ライブラリ化されたものもあります。
量子機械学習　には、古典機械学習の一部の計算を量子アルゴリズムで高速化するものと、学習モデルそのものを量子化する（量子ゲートをニューラルネットワークとみなす）ものがあります。
どちらのことを言っているのかが一見してわかりにくいです。。

  • 5qubitマシンで実行したとあるが、このGroverアルゴリズムは何bit使っている？

→ 2qubit。5qubitのうち2qubitのみ使用したことになる。

  • 無料で使える？

→ 無料。16qubitぐらいまでは使える。お金を払えば、最新の50qubit級が使えると思う。
50qubitでも、2^50の古典情報が表現可能なので、Peta bit 級が埋め込める。 (2^50 = (2^10)^5 = 1024^5 ~ 10^15 = 1Peta)
でも、問題を解くためにゲート操作をしていくと正解率がどんどん下がっていくので、答えが得られるかどうかはまだまだ注意が必要。

*1:ここでは振幅の規格化はしていません。見やすさのためです。

*2:RfゲートもNOTやXORに相当する基本ゲートを組み合わせて実現しますが、どの程度のゲート規模で実現可能かは重要です。 もしRfゲートが複雑すぎると、実用上はほとんど意味がなくなってしまいます。

*3:ここでいう計算量は、オラクル、すなわちRfゲートを何回使う必要があるかを意味します。乗算回数などで定義される通常の計算量とは異なります。

★Qiskit★ Google colab上で circuit.draw が動かない話★量子プログラミング★

Qiskit

こんにちは。Kumaです。
IBMの量子プログラミング言語　Qiskit　を使い始めました。
qiskitはpythonベースで、Google colab 上で動きます。環境構築しなくていいので初心者には大変魅力的です。
しかし、Google colab では一部の機能が動作しないか、工夫が必要です。
今回は量子回路を可視化してくれる機能　circuit.draw　が動かない問題について、説明します。

circuit.draw

circuit.draw は、量子回路 circuit に対して .draw をすることで可視化してくれる機能です。
print(circuit)でもアスキーアートが見れますが、Google colab 上では位置ずれするのであまり使えません。

しかしcolab上で、ふつうに実行するとエラーが出ます。同じ質問がコミュニティにもあります。
Google Colab - ImportError: The class MatplotlibDrawer needs pylatexenc
quantumcomputing.stackexchange.com

原因は？

原因は、platexenc がインストールされていないことです。
pip install platexenc をすればよいです。
しかし、circuit.drawを実行する前にやっておかなければなりません。
もし先に.drawを一度でもしてしまうと、仮想マシンの再起動をしないと動きません。
Colab上で　ランタイム　→　ランタイムを出荷時にリセット　をしてください。
その上で、必ず platexenc のインストールを先にしてから.drawを実行です。
ここが最初にハマったポイントでした。

.draw('latex')ではだめ、.draw('mpl')をしよう

platexencがインストールされていても、.draw('latex)は動きません。
これは別途ツールのインストールが必要なためです。Colab上で入れる方法が私には分かりません。
でも、.draw('mpl')のほうは動きます！　これはmatplotlibさえインストールしておけばOKです。
以下はqiskitのQAOAのチュートリアル　https://qiskit.org/textbook/ch-applications/qaoa.html　で回路可視化したものです。
*1

f:id:phymath1991:20201029173344p:plain

結論

まず、pip install platexenc　次に　.draw('mpl')　です。
環境構築のハードルが減っていけば、私のようなド素人でもプログラミング試しやすくなるので、期待したいです。

*1:ちなみに、QAOAチュートリアルで実機に投げるところが　ibmq_essex となっていますが、このマシンは20201029現在は存在しません。 性能向上に伴ってマシンが入れ替えられたものと思われます。例えば　 ibmq_santiago などに投げてやってください。

KumaのLTまとめ "数値積分を1000倍高速にした話"　(第18回　日曜数学会 )

こんにちは。Kumaです。
第18回の日曜数学会で発表した
「数値積分を1000倍高速にした話」
について、内容の振り返りとQ/Aで返しきれなかったものをまとめます。

Title

f:id:phymath1991:20200628235450p:plain

Agenda

f:id:phymath1991:20200628235455p:plain

モチベーション

f:id:phymath1991:20200628235459p:plain

補足説明*1

なぜ計算時間がかかるか？

f:id:phymath1991:20200628235511p:plain

モンテカルロ法による高速積分

f:id:phymath1991:20200628235518p:plain

f:id:phymath1991:20200628235523p:plain

f:id:phymath1991:20200628235528p:plain

f:id:phymath1991:20200628235536p:plain

量子アルゴリズムによる二次加速

f:id:phymath1991:20200628235541p:plain
補足説明*2

f:id:phymath1991:20200628235544p:plain
補足説明*3

f:id:phymath1991:20200628235551p:plain
補足説明*4
補足説明*5

まとめ

f:id:phymath1991:20200628235555p:plain

f:id:phymath1991:20200628235559p:plain

reference

量子コンピュータを用いた数値積分
https://www.google.com/url?sa=t&rct=j&q=&esrc=s&source=web&cd=&ved=2ahUKEwiB--zd46TqAhUCO3AKHTdZAiAQFjAAegQIAxAB&url=http%3A%2F%2Fwww.unisys.co.jp%2Ftec_info%2Ftr90%2F9003.pdf&usg=AOvVaw3zBnW66V3Ea2Crs80fExGf
Qunatum Native Dojo
Welcome to Quantum Native Dojo! — Quantum Native Dojo ドキュメント

Q&A

ニコ生コメントへの回答です

  • q#使ってます？

普段はMatlab使いなので、Matlabで書いています。
ライブラリなど無いので、自作しています。行列計算ができれば出来ます。（指数的にメモリを食いますが）
最近はpythonでIBMのqiskitも使い始めました。

  • どこまでの積分精度を必要としましたか？

相対誤差で1e-4程度でしょうか。これ以上追い求めても通信路の設計には不要だったので。

もし被積分関数が特定の構造を持てば、圧倒的に早い方法があります。
ガウス求積法が最もポピュラーです。
ガウス求積法は、通信では相互情報量の計算に使っています。どういうことかというと・・・
通信路で最も基本的なモデルとして、加法的白色ガウス通信路(AWGN)があります。
こいつの”尤度”はp(y|x) = exp(-(x-y)^2)の形をしています。
相互情報量の計算は、この”尤度”を含む積分になります。
さて、ガウス求積法によればexp(-x^2)を被積分関数に含むものは
高速に実行可能です。（わずか数点の直交多項式のy座標を求めれば良い）
ゆえに、相互情報量の計算はガウス求積法が適応可能でモンテカルロよりずっと早く計算できます。
（実装して確認しています。というか実務でバリバリ使っています）
もちろん、次元があまりに高いAWGNではモンテカルロのほうが早いかもしれません。

量子アルゴリズムが二次加速であるという事実は次元に依存しません。
つまり何次元でも量子アルゴリズムのほうが”早い”です。

  • 準乱数法*6は試しましたか？

試していませんが、もっと高速化したいので興味は持っています。

*1:実は光ファイバー中の分極現象は、通信路容量の上限を支配します。通信路容量を上げたければ、 信号入力パワーを大きくし、雑音に負けないようにすればよいです。しかしファイバー中では、分極現象が信号入力パワーが大きいほど強く発生します。 この分極たちは余計な電磁波を放射し、信号入力である電磁波にとっては雑音にしかなりません。結局、通信路容量の上限は分極の強さだけで決まります。

*2:まぁ冷却器も広い意味で量子コンピュータの一部といえますが、しかし計算の心臓部ではありません。

*3:かなり雑な説明をしています。コンセプトは合っているが、計算式としては嘘です。Referenceを参照のこと。

*4:振幅増幅はGroverアルゴリズムの一般化といえるものです。ここではQunatum　Native Dojo からGroverアルゴリズムの説明を引っ張ります

*5:振幅増幅を用いて、確率を低計算量で推定できる方法は以下のようなものです。まず確率pの小数点第1位がわかったとします。 次に小数点第2位を推定します。このとき、誤差が1/10にならないといけないわけですから、ふつうにやれば10^2倍の測定回数が必要です。 確率pの小数点第1位はわかっているので、これをもとの確率pから差し引き、これを推定する問題と考えます。この問題にたいして再び前スライドのような”都合の良い”量子状態を構成します。 そして確率増幅アルゴリズムにより|00>の確率pを10倍します。 p-> 10*p 。この状況で小数点第1位を求めます。求まった値は、元の確率pでいうところの小数点第二位です。 しかし必要な測定回数は、はじめに小数点第1位を求めたときとおなじになっています。 この手続を繰り返すことで、計算量を増加させずに任意の小数点第n桁が得られます。この方法で、必要な測定回数の総和を計算すると、なんと二次ではなく一次で済んでいます。 計算の詳細はReferenceを見て下さい。

*6:あえて偏りのある乱数を使ってモンテカルロすることで、計算の収束が早くなることがある　という話だったと思います

スマブラ物理部　ふっとび軌跡 (第15回　日曜数学会 黄黒氏リスペクト)

こんにちは。Kumaです。
第15回の日曜数学会議で黄黒氏が発表した
「大乱闘スマッシュブラザーズにおけるふっとび軌跡は空気抵抗が考慮されている」
をふまえて、空気抵抗係数に対して軌跡がどのように変わっていくのか？　という疑問を持ち
計算してみました。

背景

超有名ゲーム　大乱闘スマッシュブラザーズ　の世界では、ふっとび軌跡は放物線となりません。
ある水平位置で急速に落下するような軌跡を取ります。
この軌跡は、（パッ見でも）空気抵抗がある場合の軌跡によく一致しています。
どうも空気抵抗が実装されているらしい。

これが　第15回　日曜数学会 で黄黒氏が発見した事実です。

空気抵抗がある場合とない場合

計算結果を示します。
なお計算ツールはMatlabで、
空気抵抗が無い場合：
x = m*v0/k*(1-exp(-k*t/m))*cos(theta);
y = (m/k)*( (v0+sin(theta)+m*g/k)*(1-exp(-k*t/m))-g*t)+h0;
空気抵抗がある場合：
x0 = v0*t*cos(theta);
y0 = (-1/2)*g*t.^2+v0*t*sin(theta)+h0;
となります。
kは空気抵抗係数です。

f:id:phymath1991:20190706130816p:plain
空気抵抗がある場合(青)とない場合(黒)

原点付近では、空気抵抗があると、無いときよりもy座標が上に来るというのは面白いですね。
また、空気抵抗があるとあるx座標より先に行けずに急速に落下します。(指数的に落ちます）
この位置を測ることでスマブラ物理のパラメータを楽に同定できることも、黄黒氏は指摘しています。素晴らしい！

なお、空気抵抗がある場合の式は k = 0 を代入すると空気抵抗が無い場合の式に一致・・・しません！
k = 0 を入れると分母が0となって発散してしまいますね。
k = 0 付近では変な挙動をしているのだろうと予想できます。（真面目に近似式を得たいならたぶんTaylor展開です）

空気抵抗係数を変えた時の挙動

空気抵抗係数kをだんだんと0に近づけていってみました。

f:id:phymath1991:20190706132833p:plain
空気抵抗係数を変えた時の挙動(v_init = 2000 m/s)


以下は考察です。

kを0に頑張って持っていっても、どうやら空気抵抗が無い場合の形状にはなかなか一致しない（収束が遅い）ことが分かります。
k = 1e-10 のときの軌跡をみてみると、ジグザグしています。これがツールで扱える有効桁数を超えてしまっていることを示しています。
（切り捨て＝四捨五入が発生し、ジグザグしているわけです）
もう少し収束が早い式でないと、シミュレーションでみるには辛いわけです。

これを裏付けるために、初速度を1/10倍にしてみます。

f:id:phymath1991:20190706133024p:plain
空気抵抗係数を変えた時の挙動　(v_init = 200 m/s)

今度は、収束した様子がみられますね。
係数の調整により、だいぶ収束が早まったということですかね。


計算してみると、思っていたより大変な点が見つかって面白かったです。

まとめ

空気抵抗がある場合の曲線から、k→0の極限をとると、空気抵抗が無い場合の曲線に収束するが、
初期値によってはなかなか収束の様子が見にくい。。


それでは。