そこ曲がったら、量子坂？(上り) part4
極めて曖昧な概念でありながら、量子コンピューティングのレゾン･デートルに祀り上げられてしまった量子優位性界隈が、2025年になって喧しくなったのは、当然の帰結と思われる。誤り耐性量子コンピューター(FTQC)でも、超2次の量子加速を示す実用的なアプリケーションを見つけることは困難である・・・という予想は、量子コミュニティで広く共有されていながら、不都合な真実として、大衆向けには長らく伏せられていた(!)。それが誇大宣伝を招く遠因にもなっていたが、量子冬の時代が到来することを避けるためなのか、2025年は、｢超2次加速を示す実用的アプリケーションを見つけることは困難｣であることを堂々と主張し始めた。そして、実用的アプリケーションであれば、何らかの意味で古典を超えたら、それで良いのではないか?というキャンペーンも、大きく展開されるに至った。
何らかの意味で古典を超えるという低いハードルであれば、FTQCを待たずとも、NISQでも可能かもしれない・・・ということで、｢NISQでも役に立ちますよ｣キャンペーンは、24年後半あたりから展開されてきた。実用的アプリケーションにおいて、NISQが何らかの意味で古典を上回るという主張を行った論文。本稿で紹介する論文も、その文脈で位置付けることができる。
仏Qubit Pharmaceuticals🐾1他🐾2の研究者は、｢コンピュータ支援創薬(CADD)において、NISQ実機の実用性を確認した｣と主張する論文(以下、本論文[*241])を発表した(25年12月11日@arXiv)。しかし、この主張は危うい、と思われる。
🐾1 量子S/Wスタートアップ。
🐾2 豪Q-CTRL(S/Wスタートアップ)、仏ソルボンヌ大学。
【1】本論文の主張
本論文は、以下のように主張する：
『3次元タンパク質の水和部位予測🐾3において、NISQ実機🐾4は、量子有用性❚補足1❚を示すことが期待できる(☞【3】(3))』。
🐾3 当該予測タスクは、実用的で商業的にインパクトがあるタスクであると、本論文は考えている(☞【5】(1))。
🐾4 IBMのNISQマシン各種(つまり量子ビット･モダリティは超伝導)。量子ビット数は最大で156。
❚補足1❚
計算の高速化に限定せず(というより、高速化以外に目を向け)、何らかの意味で、量子計算が古典計算よりも優れている場合を、量子有用性と呼んでいる。先走って言うと、以下のことも本論文で示されている：本論文の量子的手法と同等若しくは同等以上に優れたな、古典ヒューリスティクス(☞【4】(1))あるいは、機械学習ベースの手法(☞【4】(4)2⃣)が存在する。つまり、主張は危うい、と思われる。
【2】事前整理
(1) タンパク質とリガンドの相互作用における水分子の配置予測(水和部位予測)
1⃣ 背景：そもそも論
本論文は、(量子)コンピューター支援創薬を扱っている。創薬モダリティは低分子化合物である。創薬プロセスでいうと、リード化合物🐾5のスクリーニングや、リード化合物の最適化(リード最適化)を効率化する(時間を短縮する)ことを目指している。別の見方をすると、生体内の標的タンパク質と低分子化合物(リード化合物)とのドッキング･シミュレーションに焦点があてられている。リード化合物は、タンパク質との結合という文脈では、リガンドということになる。
生体内のタンパク質は、親水性のアミノ酸が表面に現れるように(自発的に立体構造を)形成するので、水に易溶らしい。
タンパク質の熱力学的諸性質は周辺水分量に依存して変化するので、水の存在は、その生理学的機能発現に必要不可欠である。加えて水は、タンパク質とモダリティが低分子であるリード化合物との相互作用にも大きな影響を与える。タンパク質とリガンドの相互作用(結合)を正確に推定するには、水分子の配置(分布)を予測する必要がある。
🐾5 選抜オーディションに勝ち抜いたヒット化合物をリード化合物と呼ぶ。モダリティは低分子化合物ということになる。ヒット化合物は、薬理活性が確認できただけの化合物である。オーディションでは、医薬品として活躍できる素質(毒性がないこと、低分子化合物＝経口薬なので溶解度が高いこと、薬物動態に優れていること)が問われる。
2⃣ 従来の取り組み
タンパク質の水和構造予測手法には、実験的手法と計算的手法の二つに大別される。前者はさて措き🐾6、後者には、分子動力学(MD)シミュレーション、3D-RISM(☞下記(2))、経験的水和分布関数(empirical hydration distribution function：EHDF)を使う方法、機械学習ベースの手法(主に、畳み込みニューラルネットワーク)等がある。
MDシミュレーションは、複雑な生物系における水-タンパク質及び、水-リガンド相互作用を正確に捉えることを可能にし、薬物結合自由エネルギーの予測的な推定を可能にする。しかし、これには多大な計算コストが伴い、最先端のシミュレーションは、大規模並列処理とGPUアクセラレーションに大きく依存している🐾7。
🐾6 ㊀低温X線結晶構造解析、㊁クライオ電子顕微鏡を使った解析、㊂核磁気共鳴装置を使った解析等がある。㊁は㊀に比べて劣る、とされる。
🐾7 MDシミュレーションを量子的に高速化するタスクは、誤り耐性量子コンピューターに期待されているタスクである。
3⃣ そして量子(あくまでNISQ)への期待
NISQマシンでMDシミュレーションを高速化することは期待できない。代わりに、統計力学に基づく積分方程式論である3D-RISMに対して、NISQ実機ができることを提案した論文が、本論文の先行論文[*242]である🐾8。[*242]は、3次元の標的タンパク質を2次元にスライスして、3D-RISMの分布関数を2次元的に変換していた。このため、3次元への適用には課題があった。本論文では、3次元対応を可能としている。
🐾8 筆頭著者は、仏Pasqalの研究者。Pasqalは、量子H/Wスタートアップで、モダリティは中性原子。
(2) 3次元-相互作用点モデル(3D-RISM)の一般論[*243],[*244]
0⃣ 液体に関する基礎知識を整理
㈠ 構造
液体には流動性があるため、何らかの”構造”をイメージすることは難しい。大域的には全く、その通りである。ただし、原子間距離程度のスケールでは、原子間の相互作用により、何らかの構造(≃秩序🐾9)が現れることが期待できる。個々の原子には自由度があるので、自由度を潰す｢何らかの平均操作🐾10(＝数学的には何らかの積分計算🐾11)｣を行って初めて、(原子からなる)分子の集団に構造が定義できる。
以上をまとめると、次のようになるだろう：統計力学🐾12を活用すると、液体に構造が定義できそうである。具体的には、分布関数を用いることで、液体に構造を定義することができる。分布関数とは、密度場に対する｢モーメント｣を指す。確率･統計の文脈では、1次モーメント＝平均、2次モーメント＝分散、･･･であった。液体構造という文脈では、 1次モーメント＝密度場の平均、2次モーメント＝異なる2つの空間位置における｢密度場×密度場｣の平均、である。 2次モーメントは、2体密度分布関数🐾13と呼ばれる。2次モーメントを、1次モーメントの2乗で”正規化”した関数は、2体相関関数と呼ばれる。
🐾9 当然、短距離秩序ということになる。
🐾10 具体的には、統計力学の文脈における｢アンサンブル平均｣をとる。このように、統計力学的に由緒正しき御作法で平均操作を行うという枠組みにより、液体構造→分布関数(≃2体相関関数)から、熱力学量を自然に導出できる(ので、物理的には美しい)。
🐾11 そのため、数学的には、積分方程式理論と呼ばれる。
🐾12 多数のメンバーからなる集団に対して、そのような平均操作を行った上で、物理を展開する学問分野は、統計力学であった。
🐾13 言葉を替えると以下の通り：2体密度分布関数ρ(2)の引数をAとBとすると、ρ(2)は、A及びBにおいて分子(の中心)が同時に存在する、相対確率を表している。
㈡ 分類[*245]
物質を構成する原子間の相互作用(結合)で、液体を分類することができる。この場合液体は、㊀分子性液体、㊁イオン液体、㊂金属液体、の３つに分類される。㊀分子性液体間の結合は、水素結合あるいはファンデルワールス力などの分子間力であり、弱い。㊁イオン液体間の結合は、陽イオンと陰イオンの間に働くクーロン力で、比較的強い。㊂金属液体間の結合は、金属結合であり、強い。
水やアルコール･有機溶媒といった”普通の”液体は、分子性液体である。
1⃣ RISM
0⃣㈠において展開した、液体構造を定義する枠組みは、数学的には、積分方程式理論と呼ばれるのであった。積分方程式理論を、多原子分子で構成される分子性液体に適した形式に調整したモデルが、相互作用点モデル(Reference Interaction Site Model：RISM)である🐾14。
分子性液体を対象とする場合、分布関数(≃相関関数)を記述するためには、分子の位置と配向を考慮しなければならない。RISMでは、溶質分子と溶媒分子の両方の｢配向｣を平均化して自由度を潰してしまう。このためRISMでは、溶質分子の位置と溶媒分子の位置のみが自由度として残り、分布関数(2体相関関数)は、距離だけの関数(つまりは、動径分布関数)になる。ここで、距離rは、r＝｜溶媒分子の位置ー溶質分子の位置｜である(｜･｜はL1ノルム)。液体が均一系であれば、記述するモデルとしては、RISMで十分ということになる。
RISMはミニマルである(＝情報を縮約し過ぎている)ため、タンパク質とリガンドの結合自由エネルギーを求めることはできない(という理解)。つまり、創薬支援という文脈では、RISMは力不足と言えるだろう。
🐾14 不思議なことに、Referenceには特に意味がないらしい。このため、和訳には現れないことも多い。つまり、"参照”相互サイトモデル、という和訳でなくとも良いらしい。
2⃣ 3D-RISMについて概要→詳細は、(3)へ(⤵)
分子科学研究所の平田文男名誉教授他によって開発されたモデル。3D-RISMでは、溶質の配向は残し、｢溶媒の配向｣のみを平均化して潰す。RISMでは、分布関数gはスカラー量である距離rの関数であった：g=g(r)。3D-RISMでは、rが溶質分子の配向自由度を反映し、ベクトルrになる。つまり、g=g(r)となる。分布関数は溶媒分子の原子毎に得られる。つまり、水の場合、水素原子と酸素原子の分布関数が得られる。通常は、酸素原子の分布関数(のみ)を水分子の分布関数と見做すようである。
(3) 創薬支援という文脈における3D-RISM
1⃣ メイン
(2)2⃣の分布関数から、熱力学量が自然に導出される🐾15。創薬支援の文脈では、結合自由エネルギーが、重要である🐾16。より正確に言えば3D-RISMでは、2つの溶質(標的タンパク質とリガンド＝モダリティが低分子であるリード化合物)の脱溶媒和自由エネルギーを、分布関数から算出することができるので、創薬支援に資する、ということができる。
リード化合物のスクリーニングあるいはリード化合物の最適化(リード最適化)には、リード化合物とターゲット(＝タンパク質)との結合自由エネルギーを(できるだけ正確に)予測する必要がある。3D-RISMでは、タンパク質とリード化合物(リガンド)の分布関数を求めることで、結合自由エネルギー(脱溶媒和自由エネルギー)を求めることができる。
創薬と言う文脈＆タンパク質と(リード)化合物というセットアップで、結合自由エネルギーを算出する場合、厳密に言えば、水中のタンパク質の｢構造揺らぎ｣を考慮しなければ、正確に算出することはできない。それは、タンパク質が最安定構造の周りで時間的・空間的に揺らいでいるから、と説明される。残念ながら3D-RISMは、溶媒(水)における溶質(標的タンパク質)の｢構造揺らぎ｣を考慮できない。従って、 構造揺らぎをなんらかの方法(例えば、分子動力学法)で補った3D-RISMが創薬支援では、使用されている(という理解)。
🐾15 自然に導出されるものの、あくまで3D-RISMは近似理論なので、分布関数も近似的にしか求めることができない。従って、熱力学量も近似的しか求められない(ので、当然、誤差が含まれる)。また、3D-RISM により得られる分布関数は、結晶構造解析から得られる水和構造を再現しないことが指摘されているらしい[*246](☞【5】(2))。
🐾16 結合自由エネルギーは概ね、㊀2 つの分子の結合前後の原子間相互作用エネルギー変化、㊁溶質分子のエントロピー変化、㊂化合物(リガンド)とターゲット分子(標的タンパク質)の脱溶媒和自由エネルギー、の3 つの寄与に分けられる[*247](らしい)。
2⃣ 溶質分子周りの溶媒分子の分布関数を求めるプロセス
数学的には、RISM及び3D-RISMは積分方程式論であった。つまり、3D-RISMで溶媒分子の分布関数を求めるには積分方程式を解く必要がある。水和させる溶質(標的タンパク質)の分子構造における各原子の3次元座標🐾17が入力として必要となる。分子構造(3次元座標)に加えて、
3D-RISMは力場パラメータも必要とする。本論文では、ff14SB Amber力場🐾18を使用し、溶媒(水)の記述にはSPC/E 水モデル🐾19を使用している。
3D-RISMシミュレーションは標準CPUで実行され、平均シミュレーション時間は約30分だったらしい。
🐾17 PDB(Protein Data Bank ：タンパク質構造データバンク)から取得している。PDBは、タンパク質、核酸、糖鎖など生体高分子の3次元構造の原子座標(立体配座)を蓄積している国際的な公共のデータベースである。
🐾18 複雑なことはさて措き、Amber(☞🐾20)では最も代表的な力場の一つらしい。水モデルとして、ff14SBはTIP3P(⤵🐾19)を使用している。本論文の水モデルはSPC/Eであるが、整合性が付くのか、よく分からない。出典：https://www.ag.kagawa-u.ac.jp/charlesy/memo/md-simulation/
🐾19 代表的な水の分子動力学シミュレーション･モデルには、SPC/E(Extended Simple Point Charge)とTIP3P(Transferable Intermolecular Potential 3-point)があり、生体分子シミュレーションではTIP3Pが広く使用されている(らしい)。SPC/EとTIP3Pは(水分子を構成する)酸素原子と水素原子に、1つずつ｢相互作用点｣を配置して、3点からなる剛体として扱う(らしい)。出典：https://www.tus.ac.jp/today/archive/20230914_4819.html
【3】本論文の技術要素
(1) 水和部位予測のワークフロー
1⃣ メイン
本論文は、量子コンピューターを用いたタンパク質の水和部位予測をタスクとしている。当該タスクのワークフローは、以下の通りである。まず、3D-RISMを使って、標的タンパク質(溶質)周りの水分子(溶媒)の分布関数g(r)を取得する🐾20。なお、この文脈での水分子は、水和水分子と呼ばれるし、そう呼ぶほうがしっくりくる。次に分布関数g(r)を、｢正規分布の線形結合｣として表現する。ここで、正規分布の平均は、水和水分子の位置である。分散は、位置の不確実性である。
なぜ、わざわざ、｢g(r)を正規分布の線形結合として、表現するのか？｣は、(2)で触れる。いずれにしても、本論文では、こうして『水和水分子の位置を推定する問題』＝『水和部位予測問題』を、g(r)を近似する”最適な”混合正規分布を見つける｢最適化問題｣に変換する。最後に、この｢最適化問題｣をQUBO(2次制約なし二値最適化)問題に変換する。QUBO問題に変換することで、量子最適化のフレームワークに持っていき、量子有用性の期待を抱かせている。
🐾20 g(r)の計算には、無料で使用できるAmbertoolsパッケージ(の3D-RISMツール)を使用している。Ambertoolsは、生体分子シミュレーション･ソフトウェアであるAmberの関連ツール群。3D-RISMも含まれている。なお、g(r)のrは3次元座標を表すベクトル。
2⃣ 追記
本論文は、あくまでゲート方式の量子コンピューターを想定している。つまり、QUBOを量子アニーリング･マシンで解くわけではない。IBMのNISQ実機(モダリティは超伝導)上で豪Q-CTRL(S/Wスタートアップ)の量子最適化ソルバーを用いて🐾21、QUBO問題を解いている🐾22。量子アルゴリズムは、QAOA(近似的量子最適化アルゴリズム)ということになる。
🐾21 Q-CTRLのパフォーマンス管理ソフトウェアFire Opalを使った量子誤り抑制を実施している。Fire Opalはパルス制御等を通じて、量子プロセッサのハードウェア･ノイズやデコヒーレンスを極小化することで、量子誤りを"抑制"する。
🐾22 エネルギーが最小となるような水和水分子の位置が求められる。
(2) 分布関数を混合正規分布で表現仕直すことの意味
3D-RISMだけでも、水和構造予測ができるのに、なぜ本論文では、3D-RISMで得られる分布関数を、正規分布の線形結合として表現しているのか。なぜ、こんな面倒なことをやっているのか？というと、主に２つの理由があると思われる。一つ目は、3D-RISMの制限対応と考えられる。正規分布の線形結合という算術的な表現を物理的に解釈すると、｢各水和水分子が、タンパク質周囲における特定の安定位置周りを、(正規分布に従う乱数で表現される乱雑さという意味で)ランダムに振動している｣、と仮定していることに相当する。この枠組みを採用すると、3D-RISMが(標的)タンパク質の構造揺らぎを考慮できない、という欠点を一部回避することが可能と考えられる。
2点目は、効率的な離散化対応と考えられる。3D-RISMで得られる分布関数は連続関数であり、水和部位を予測するには、離散化しなければならない。予測精度を上げるためには、細かい離散化が要請されるが、それは計算コストを上げることにも繋がる。必然的に、効率的な離散化が望まれる。本論文は、混合正規分布で離散化することで、効率的な離散化を実現していると考えられる。
(3) 量子有用性の期待に至る道筋：あらすじ
1⃣ まず、古典的手法で扱えるサイズの問題を対象とする。その対象に対して、本論文の量子最適化手法が、古典的手法と”同等の”成功確率で、最適解を見つけられることを示した(☞【4】(1))。
2⃣ 1⃣で本論文の量子最適化手法がワークすることを示した上で、問題サイズが大きくなると、CPLEXでは手に負えなくなることを、天下り的に主張する(☞【4】(2))。加えて、量子最適化手法では、問題サイズが大きくなるにつれてスケールする、つまり性能が上がることを示した(☞【4】(3)2⃣及び4⃣)。
3⃣ さらに、水和部位予測タスクをQUBO問題に変換して解くというQUBOアプローチを、他アプローチと比較して、次のように主張する：QUBOアプローチは、他アプローチと同等あるいは場合によってはそれを上回る性能を示す(☞【4】(4)2⃣)。
4⃣ まとめると、QUBOアプローチは他アプローチよりも予測性能が高い上に、スケールする。従って問題が大規模化すると、量子コンピューター上で実行することにより、量子有用性が現れることが期待できる。ただし、シミュレーテッド･アニーリングのような古典ヒューリスティックや機械学習ベースの手法より、予測性能が高いか否かは、未解決である(という理解)。
【4】評価結果、比較結果
(0) 評価指標
水和部位予測を目的とした数値手法の品質を評価するために、結晶水(CW)🐾23の位置を使用するのが一般的である。正確なCWの位置(正解)は、実験データで利用可能である。
具体的には、4つの評価指標が本論文で使われている。❶C、❷⟨CS⟩、❸P∗、➍⟨P⟩、である。❶と❷は、対象(今の場合、標的タンパク質)全体にわたり、予測が成功したか?を計測する指標である。成功とは、設定したしきい値の範囲内に、正解と予測値が含まれるか?で判断する。付け加えると、❶がメインで❷は補助的な指標である。❸と➍は、予測値が正解にどれだけ近いか?を計測する指標である。❸は、最良(最近接)値の算出を目的としている。➍平均的な値の算出を目的としている。細かい説明は、下記❚補足2❚を参照。
🐾23 結晶の中に、一定の比で化合物に結合している水のこと。
❚補足2❚
❶は正解の周囲3 Å 半径内に｢少なくとも 1 つ｣の予測値が存在するかどうかを探索して計算される。しきい値が3 Åということになる。誤差3 Åで予測できていれば、正解に対して1のフラグを立てる。予測できていなければ、0にフラグをたてる。全ての正解に対する、｢１のフラグがついた正解｣の割合という理解で良い。
❷は、誤差3 Åで複数の予測値があれば、その個数をそのままカウントする。❶は、複数の予測値が見つかっても、１のフラグを立てるのみであった。❷＝(正解と｢しきい値範囲｣内にある予測値の個数)/(フラグ1が付いている正解の総個数)、である。従って、❷は1を超える可能性もある。あくまで補助的存在である。
❸は、予測値と正解との間の、最も近い距離を表象する。➍は予測値と正解との距離を、標的タンパク質全体にわたって平均した値を表象する。具体的に数式もどきで表すと、❸は、∑max[1/(1+r)]という式で計算される。ここで、rは予測値と正解との距離である。同じ位置に存在した場合に０(つまり、supが0)になり、1/(1+r)＝1となる。つまり、❸でとっている最大値(supは1)は、予測値と正解とが、最も近いときに発生する。なお、∑は正解の個数にわたるので、実際は、正解の個数で割って”正規化”する。従って、1を超えない。➍は、1/(1+r)の平均値を計算する。ここでも正規化して、➍が1を超えないようにしている🐾24。
🐾24 本論文における➍の定義式である式(5)には、誤植があると思われる。
(1) 古典手法と同程度の成功確率でグランドトルゥースを得られる
PDBから選択した6つのタンパク質🐾25について、QUBO問題に変換した水和部位予測タスクを、3つの方法で解いた：①本論文の量子最適化手法、②シミュレーテッド･アニーリング(SA)法、③貪欲法。③は2つのタンパク質に対してのみ、最適解を見つけられたに過ぎない。加えて、グランドトルゥースとの一致する割合(成功確率)も低かった。グランドトルゥースは、CAPLEXを使って得られた解としている。
①と②の勝負は、①の4勝5敗2分け(5つのタンパク質が2パターン＋１つのタンパク質で、都合11ケース☞🐾25)と考えられる。圧勝が3ケース(3/4)、完敗が1ケース(1/5)。SAと同等と言って良いのだろう。なお、成功確率は、①で40%弱(②で45%程度)と低い。
🐾25 タンパク質の結合部位内で、L：リガンドの静電的効果と立体的効果を考慮に入れた場合と、NL：考慮に入れない場合、の2パターンを考えている。6つ目のタンパク質は、リガンドを考慮していない(アポ構造)。
(2) CPLEXの限界と量子有用性が実現する時期
QUBOベースの水和部位予測手法は、1,000変数程度の規模で実用的に意味を持つらしい。ザックリ言って、900変数が、0.5Åのグリッド間隔に相当する。正確な古典的ソルバーとするCPLEXは123変数でQUBO問題を解くのに苦労しているらしい。このため900～1,000変数で、｢古典的には解けない、実用的に意味を持つ問題が、量子では解ける｣と言う意味での量子有用性が実現する、と主張している。もっとも、ヒューリスティクスであるSA法は、900～1,000変数でも、苦労せずに解くことはできる。
併せて量子リソース推定も行っている。量子ビット数と2量子ビットゲートの数をプロットしてフィッティングした図において、量子ビット数(変数の数)を156から900 まで伸ばして、2量子ビットゲート数を推定した。その結果、約100,000個の 2 量子ビット ゲートが必要になるという結果がでた。H/W開発各社の物理量子ビット数のロードマップと照らし合わせて、この数の2 量子ビット ゲートが実現する時期は、2028～2029年と推定している。
(3) QUBOアプローチがスケールする
0⃣ 為念：改めて言葉の整理
QUBOアプローチとは、水和部位予測タスクをQUBO形式の最適化問題に変換して解くことを指している。QUBOアプローチという文言は、CPLEXで解くことや、機械学習ベースの手法で解くこととの対比で用いられる。従って、QUBOアプローチは、本論文の量子最適化手法で解くことに限定されず、シミュレーテッド･アニーリングで解く場合も含まれる。
1⃣ 概要
QUBOアプローチが123変数まで、スケールすることをNISQ実機で(つまり、本論文の量子最適化手法で)確認している。スケールするとは、変数が{44→71→123}と増加する(離散化のグリッド･サイズが小さくなる)につれて、水和部位予測の予測性能が向上することを指している。予測性能の比較には、(0)の評価指標❶～➍を使っている。ちなみに、123変数は、グリッド間隔が1.9Åである。
なお、使用したNISQ実機の量子ビット数は最大156であり、 各バイナリ変数(水分子が空間位置に存在していれば1、存在していなければ0)は、1量子ビットにマッピングされる。
2⃣ 追記１：NISQ実機の結果を定量的に言うと
変数{44→71→123}において❸≒➍であるから以下、❸と➍のセットでは、➍のみ言及する。1⃣で示した通り、 変数が{44→71→123}と増加するにつれて、❶、➍は値が大きくなる(予測性能が上がる)。❶＝{0.45程度→0.5程度→0.6程度}、➍＝{0.2程度→0.25程度→0.3程度}である。
ちなみに、❷＝{1程度→1程度→1.1程度}である。
3⃣ SAを援用して
さらに、123変数を超えて(当然ではあるが･･･)QUBOアプローチがスケールすることを確認している。ただし、156変数を超えるとNISQ実機では対応できないため、古典手法であるシミュレーテッド･アニーリングで解くことでスケールすること確認している。変数は{900→1450→6824}である。
4⃣ 追記2：SA拡張バージョンを定量的に言うと
ここでも、 900変数と1450変数では❸≒➍であるが、 6824変数では異なる。❶＝{0.82程度→0.9程度→0.9程度}。6824変数で打ち止め感はあるが2⃣と比べて、相当改善されている。➍＝{0.42程度→0.42程度→0.45程度}。❸＝{0.42程度→0.42程度→0.4程度}。➍も2⃣と比べて、改善されていることがわかる。一方で、900変数以降の改善は、ほぼない。その傾向は❸でも観察できるが、 6824変数に至っては、性能が落ちている🐾26。なお、❶と➍の計算式は異なるので、改善幅を同列に議論することはできない。ちなみに、❷＝{1程度→1.5程度→2.8程度}である。
🐾26 この件について、本論文で、考察等は特に行われていない。
(4) QUBOアプローチと他手法の比較
1⃣ 比較する他手法
⓪ QUBOアプローチでは、900変数と3974変数で、予測タスクが実行されている(それぞれ、⓪_900、⓪_3974と表記する)。つまり、シミュレーテッド･アニーリングで解いた結果である。
① Hydraprot(deep Hydration of Proteins)[*248]・・・タンパク質構造周囲における水分子の酸素原子の正確な位置を予測する。3D U-netと多層パーセプトロン(MLP)という、2つの相互接続された深層学習アーキテクチャを活用する。当該手法は、タンパク質の粗いボクセル･ベース表現を導入することから始まり、これにより3D U-netを介して候補となる水の位置を迅速にサンプリングすることができる。次に、これらの水の位置は、MLPを用いて水とタンパク質の関係をユークリッド空間に埋め込むことで評価される。最後に、事後処理ステップを適用することで、MLP予測をさらに精緻化する。
② Watgen[*249]・・・タンパク質とペプチドの原子座標が与えられた場合に、タンパク質-タンパク質若しくはタンパク質-ペプチド界面における水分子の位置を予測するアルゴリズム。
③ Placevent[*250]・・・3D-RISMから得られる連続的な分布関数から、確率密度が最大となる位置に(水分子であれば水分子の酸素原子を)配置することを繰り返して、溶媒分子やイオンの具体的な座標(離散的な位置)を予測･配置するアルゴリズム。
④ Dowser＋＋[*251]・・・タンパク質構造内の水分子の位置を予測するDowserを半経験的に改変したアルゴリズム。AutoDock Vina🐾27およびWaterDock🐾28の使用が前提となっている。
①～④は、すべて計算コストが低く、平均して通常の CPU で数分の計算しか必要としない(らしい)。☞【5】(0)
🐾27 分子ドッキングとバーチャル･スクリーニングのためのオープンソース･ソフトウェア。
🐾28 WaterDockは、タンパク質構造中の水分子の位置を予測する手法。AutoDock Vinaを用いて水分子を構造にドッキングし、スコアの低い部位を除外し、残りの部位をクラスター化する。クラスターの重心が予測される水部位となる[*252]。
2⃣ 比較結果
まず、総論をザックリ言うと、QUBOアプローチより予測性能が高い｢1⃣他手法｣が存在する。以下、個別に述べる。まずは定性的に。①は、❶,❸,➍の全てで、QUBOアプローチよりも優れている(評価指標の値が大きい＝良い)。②は、❶のみQUBOアプローチよりも優れている 。なお❶だけであれば、②が最も優れている。③及び④は❶,❸,➍の全てで、QUBOアプローチよりも劣る。❷は定量的についてのみ、言及する。
定量的には、以下の通り：⓪_900＝{❶＝0.82程度,❸＝0.42程度,➍＝0.42程度}、⓪_3974＝{❶＝0.95程度🐾29,❸＝0.5程度🐾30,➍＝0.45程度}、①＝{❶＝0.98程度,❸＝0.62程度,➍＝0.6程度}、②＝{❶＝1.0程度,❸＝0.5程度,➍＝0.4程度}、③＝{❶＝0.25程度,❸＝0.1程度,➍＝0.1程度}、④＝{❶＝0.35程度,❸＝0.2程度,➍＝0.2程度}。
なお、❷＝{⓪_900＝1.0程度、⓪_3974＝1.7程度、①＝1.2程度、②＝3.2程度、③＝1.0程度、④＝1.3程度}。
本論文は、甘めに、次にように結論している：他の古典的な手法と同等、あるいは場合によってはそれを上回る性能を示す。
🐾29 0.9を超えて、0.95程度はあるように見える。(3)4⃣と、やや矛盾するかもしれない。
🐾30 この結果も、(3)4⃣と、やや矛盾するかもしれない。
【5】考察
(0) 既述通り、3D-RISMシミュレーションの実行時間は、CPUで約30分。故に(3D-RISMの結果を使う)本論文の量子手法は(量子計算が如何に速くとも)、30分以上の時間を要する。一方で、①Hydraprotは、CPUで数分。それで、①は本論文の量子手法(⓪_900と⓪_3974)を上回っている。普通に考えると、勝負にならない。
NISQでなく後期FTQCであれば、3D-RISMを使わず、MDシミュレーションによる水和部位予測を(高い予測性能を維持したまま)高速化することが期待されている。つまり、後期FTQCであれば、機械学習ベースの手法を上回るだろうと期待されている。
(1) 本論文の創薬モダリティは、低分子化合物である。故に、実用的で商業的にインパクトがあるアプリケーションと言っていいのか？という疑問は当然ある。しかし、抗体薬物複合体(Antibody-Drug Conjugate：ADC)🐾31を設計するには、低分子化合物の知見が必要であることが知られている。具体的に言うと、結合剤(リンカー)の設計に、水和部位の予測が役に立つと考えられる。結論として、現代(2020年代)の創薬においても、水和部位の予測は有用と考えられるだろう。
🐾31 代表的なADCとしては、第一三共のブロックバスター抗がん剤であるエンハーツがある。
(2) 🐾15でも触れたが、3D-RISM により得られる分布関数は、結晶構造解析から得られる水和構造を再現しないことが指摘されているらしい[*246]。しかし、これは、 3D-RISMを3D-RISM＋MD(分子動力学)にする等の方法で対処することが可能であろう。また本論文の枠組みでは【3】(2)でも触れたように、 ｢3D-RISMでは、溶媒(水)における溶質(標的タンパク質)の”構造揺らぎ”を考慮できない｣という欠陥に対処していると考えられる。
結論として、本論文の枠組みは｢3D-RISM により得られる分布関数が、結晶構造解析から得られる水和構造を再現しない｣という問題に、対応できている(あるいは対応することは可能)と考えられるだろう。つまり、枠組みは破綻はしないだろう。
(3) 本論文手法がシミュレーテッド･アニーリング(SA)を超えるか否かは、未解決である。そもそもヒューリスティクスは、なぜ高性能なのかを理論的に説明できない手法なのだから、実際にガチンコで結果を突合して、比較するしかない。言葉を替えると、タスクによってSAを超えたり、超えなかったりするわけである。ただ、量子アニーリングとは、変数の数を(123より)もっと増やして比較することができるはず。その比較結果には、興味がある。NISQでも役に立ちますよキャンペーンでありながら、代替H/Wが存在するというのは、打ち手としてどうなのだろうか。
なお、3D-RISMの分布関数を基にしたPlacevent(③)と本論文手法を比較すると、本論文手法が圧勝するが、それをもって、NISQの実用性が確認されたということにはならないだろう。
(4) 本論文手法は、3D-RISMの分布関数を近似する混合正規分布を探索するタスクを、QUBO問題を解くというタスクに変換している。この枠組みにおける混合正規分布の探索は、ガウス過程回帰(GPR)でも同じこと出来ると思われる。Pros＆Consで両者を比較することは、意味があるように思われる。GPRはCPUでも可能だが、QPUをつかうことにより、近似性能(→予測性能)が上がる可能性があるだろう。
(0) 前置き
量子化学シミュレーションは、(ミクロ･スケール及びマクロ･スケールの)物理シミュレーションと並んで、”実用的な”量子優位性が最も期待されるアプリケーションである。ただし、量子コンピューティング✖量子化学シミュレーションの文脈において、魅力的な対象は限定的とされる[*A-77](☛詳細は、こちらを参照)🐾1。そのため、この文脈で紹介される分子は、FeMoCo🐾2やシトクロムP450🐾3等、判で押したようにワンパターンである。
そこに突然、光線力学的療法(PhotoDynamic Therapy：PDT)の光感受性物質という、期待の新人が現れて、一筋の光明が差し込んだらしい。IQM🐾4は『Wellcome Leap Q4Bio Challenge🐾5の一環として、Algorithmiq🐾6は(IQMの量子処理ユニットを用いて)、量子シミュレーションが次世代のPDTの実現にどのように貢献できるか？を実証した』と自社ブログ[*253]にてアピールした(25年7月28日)。Algorithmiqは自社サイトにおいて、｢PDTに最適な薬剤を特定する量子ブーストAIプロトコルを開発している｣と述べている[*254]。具体的には、｢能動学習フレームワークを活用して、化学空間を効率的に探索する、ハイスループット仮想スクリーニング･パイプラインを開発している｣らしい。スクリーニングの品質を維持しながら、計算コストを大幅に削減できる、と主張する。
🐾1 ちなみに、量子物理シミュレーションであれば、魅力的な対象は豊富とされている。
🐾2 農業における窒素固定を担う｢土壌生態圏における生物窒素固定系｣で中心的な役割を果たす、真正細菌･古細菌が持つ窒素固定酵素｢ニトロゲナーゼ｣の活性中心に含まれる｢鉄とモリブデンを含む補助因子｣。
🐾3 一般的に、肝臓において解毒を行う酵素として知られている。
🐾4 フィンランドの量子H/Wスタートアップ。モダリティは超伝導。
🐾5 ヒト健康分野への量子コンピューティング適用を”加速すること”、を目指しているチャレンジ･プログラム。バジェットは、US$40mil＋10mil。
🐾6 フィンランドの量子S/Wスタートアップ。創薬が主要ターゲット。
(1) 本題
それから凡そ5ヶ月後、加Xanadu🐾7の研究者は、｢誤り耐性量子アルゴリズムを用いて、PDTに有望な光感受性物質候補を特定する方法｣を開発･提示した論文(以下、本論文[*255])を発表した(25年12月17日@arXiv)。Algorithmiqは、量子古典ハイブリッド手法を採用し、NISQマシンで出来ること(今出来ること)を追及した。やや厳しい表現を採用すると、｢治療効果が十分とは言えない現状レベルに留まる｣光感受性物質のスクリーニングを高速化している、に過ぎない。
一方Xanaduは、現状レベルの光感受性物質ではPDTの効果は不十分であり(＝他の療法と比較して競争力がないので)、十分な効果を発揮する光感受性物質を探索･特定する方法が必要という立場をとる。そして、そのような光感受性物質の探索において、精度とスケーラビリティを両立させるには、量子計算(量子アルゴリズム)が不可欠であると主張する。本論文では、該当する量子アルゴリズムを実際に開発･提示し、併せて量子リソース推定を行った。リソース推定の結果、開発した量子アルゴリズムを走らせるには、誤り耐性量子コンピューター(FTQC)が必要であろう、と結論されている。
その結論は正しいだろうが、皮肉なことに、FTQCが必要なレベルであっても、量子優位性(量子有用性)が顕現するとは限らないことも示唆されている(と思われる)。
🐾7 量子H/Wスタートアップ。モダリティ:光。米NASDAQ及び加トロント証券取引所にSPAC上場することが決定(25年11月3日)。
【1】本論文の主張
本論文は、以下を主張する：
(0) 誤り耐性量子コンピューター(FTQC)を必要とする、実用的アプリケーションを見つけた(⤵(3))。
(1) PDTの治療効果を、2つの量子力学的計算変数に置換した。
(2) 十分な精度とスケーラビリティを両立して、(1)変数を計算する量子アルゴリズムを開発した。
(3) (2)量子アルゴリズムの量子リソースを、以下のように推定した：必要な論理量子ビットは、概ね180～350個。量子ゲートは、トフォリ･ゲート❚補足1❚は107～109個必要(⤴(0))。
❚補足1❚
量子ビットを採用する場合において、トフォリ･ゲートを分解しようとすると、少なくとも 6 個の CNOT ゲートが必要である。量子トリット(qutrit)であれば、3個のCNOT ゲートで分解可能である[*256]。なお、量子トリットとは、計算基底を3つ使った量子情報単位である。計算基底が2個の場合が、量子ビット(qubit)である。
【2】事前整理
(0) 量子･コンピューター支援創薬
0⃣ プレイヤー･レビュー
創薬支援を射程に収める量子S/Wスタートアップとしては、以下があげられる：日本QunaSys、米Quantum Simulation Technologies、米Polaris Quantum Biotech、独HQS Quantum Simulations、独Quantistry、仏Qubit Pharmaceuticals、フィンランドAlgorithmiq、デンマークMolecular Quantum Solutions、デンマークKvantify、英Cambridge Quantum Computing(米クオンティニュアムの一部門)、英KUANO、加Good Chemistry Company、加ProteinQure🐾8、加Menten AI、韓国Qunova Computing。
🐾8 少なくとも一時期、日本の製薬大手、第一三共とパートナーシップを結んでいた。
1⃣ 研究成果
量子･コンピューター支援創薬に関する研究成果例としては、以下があげられる：㊀古くは米インシリコ･メディシン🐾9他の研究：低分子医薬品候補物質の発見において、量子古典ハイブリッドGAN(敵対的生成ネットワーク)が、古典GANよりも優れている(☛こちら)。㊁インシリコ･メディシン他の研究：量子古典ハイブリッド生成モデルは、古典生成モデルと比べて、ドッキング･スコアを基準として、優れたヒット化合物を探索できる(☛こちら)。㊂米デロイト･コンサルティング他による研究：薬理活性)予測において、量子古典ハイブリッド･カーネル法は古典カーネル法よりも優れている(☛こちら)。㊃仏Qubit Pharmaceuticals🐾10他の研究：3次元タンパク質の水和部位予測において、NISQ実機は、量子有用性を示すことが期待できる(☛こちら)。
🐾9 米国のAI創薬スタートアップ。アルファ･フォールドを使って、肝細胞がんを対象とする新規ヒット化合物を発見。それ以前に、AIを使って発見した｢特発性肺線維症｣治療薬を開発している。
🐾10 量子コンピューティングで創薬に要する時間とコストを大幅に減らすことを目指している量子S/Wスタートアップ。腫瘍、炎症性疾患及び感染症をターゲットとしている。
2⃣ コラボレーション
㈠ QunaSys
QunaSysは、日本の中外製薬とパートナーシップを結んでいる(いた?)。新しいところでは、QuEnAIS(Quantum-Enhanced AI Synthesizer for Drug Discovery)というプロジェクトに携わっている[*257]。同プロジェクトは、2025年8月〜2027年7月に渡って行われ、デンマーク･イノベーション基金と独連邦教育研究省が資金提供する。他メンバーは独Cortex Discovery🐾11と独フラウンホーファー技術･経済数学研究所🐾12のHPC部門。
QunaSysは、テンソル･ハイパー縮約(THC)を適用したDFT計算における電子反発積分の高精度近似と、qubitization(量子ビット化)を活用した量子位相推定法(QPE)の改良を担当する。
🐾11 前臨床試験に特化した、医薬品特性予測用深層学習モデルを持つ。250万個の分子について、7億5,000万個のデータポイントを備えた、膨大な生物学実験データベースを基に開発した。
🐾12 和訳の出典はhttps://www.fraunhofer.jp/ja/institutes-establishments.html。HPC＝高性能コンピューター(≒スパコン)。
㈡ Quantum Simulation Technologies(Q Simulate)
Q Simulateは、日本にも拠点(日本支社?)がある。Q Simulate、クオンティニュアム(Cambridge Quantum Computing)及び三井物産は、量子古典ハイブリッドプラットフォームQIDO(Quantum-Integrated Discovery Orchestrator)を開発、25年8月19日から提供を開始した[*258]。
㈢ 米IBM
IBMは、米クリーブランド･クリニック🐾13と10年間におよぶDiscovery Acceleratorパートナーシップを2021年3月に締結。2023年3月には、IBM製NISQ実機｢IBM Quantum System One｣がクリニックに設置された。成果としては、IBMとクリーブランド･クリニック･ラーナー研究所が、｢タンパク質構造予測｣について研究成果を発表している(24年5月)[*129]。また、英ハートリー･センター🐾14を含めた3者は、人工知能(AI)や量子コンピューティングを含む、高度コンピューティング技術を通じて、医療と生物医学を進歩させることを目的とした共同プロジェクトを発表している(24年6月)[*130]。
🐾13 米オハイオ州第2の都市、クリーブランドにあるクリニック。Newsweekが毎年発表するよい病院の世界ランキング｢World's Best Hospitals｣の常連上位組。25年は、世界2位(1位は、やはり米国のメイヨー･クリニック)。心疾患では、全米1位と言われている。
🐾14 英国のスーパーコンピューティング･センター。
㈣ 米IonQ
新しいところでは、IonQが、ハイブリッド量子技術および量子AI技術を用いた次世代治療法開発を加速するため、カナダの再生医療商業化センター(Centre for Commercialization of Regenerative Medicine：CCRM)🐾15とのパートナーシップ締結を発表した(25年12月)[*259]。このパートナーシップは、バイオ･プロセスの最適化、疾患モデリング･ワークフロー、先進的な治療薬の設計と製造を支援する量子シミュレーション等に、焦点を当てている。
🐾15 2011年に設立されたCCRMは、加トロント大学に拠点を置くグローバルな官民連携組織。再生医療の可能性を実現し、深刻な慢性疾患に対する永続的な治療法、さらには治癒をもたらすことを目指している。特に細胞遺伝子治療に重点を置いている。加政府、オンタリオ州政府をはじめ、主要な学術機関や産業界のパートナーから出資を受けている。参考サイト：https://www.phchd.com/jp/phc/news/2025/0210
3⃣ 大規模プロジェクト
Wellcome Leap Q4Bio Challenge(https://wellcomeleap.org/q4bio/)が有名である(☛こちら)。
(1) 光線力学的療法(PhotoDynamic Therapy：PDT)
1⃣ 概要
PDTとは、光感受性物質🐾16を用いて活性酸素種(ROS)🐾17(一重項酸素)を生成し、腫瘍組織(腫瘍や新生血管)を選択的に変性壊死させる｢標的がん治療法｣である。光感受性物質は、自身の吸収帯に一致する特定の波長域の光によって活性化されるまで不活性な状態を維持する。なお、この場合の｢光｣は、レーザー光である。
PDTは、従来の他治療法よりも副作用が少ないことが一般的である。ただし、既存のPDTは、特に深部腫瘍の治療及び、高い治療効果の達成において、重要な限界に直面している。
🐾16 特定波長の光を吸収した結果、化学的若しくは物理的に活性化する物質。光を当てない限り、毒性を示さない。
🐾17 反応性の高い酸素種の総称。スーパーオキシド、過酸化水素、ヒドロキシラジカル、一重項酸素などがある。PDTでは、一重項酸素を使うらしい。
2⃣ PDTは限界に直面している
PDTの臨床的成功は、高い光感度と高いROS生成効率を持つ光感受性物質の利用可能性にかかっている。しかし、既存の多くの光感受性物質は、組織への浸透に適さない波長で吸収を示すため、体内の腫瘍を効果的に標的とする能力が制限されている。また、十分な光照射量が得られても、多くの光感受性物質は、項間交差(ISC)❚補足2❚効率が低いため、ROSの生成が制限され、治療効果が低下する。
実のところ、光感受性物質の光物理的特性は、広範囲に調整可能である。つまり、治療波長域における累積吸収とISC速度の両方を、大幅に向上させることは理論上可能である。ただし、光感受性物質の光物理的特性向上タスクには、通常、｢膨大な実験的合成、精製、及び特性評価｣作業が必要である。当該作業は、時間とリソースを大量に消費する。したがって、効率的な候補スクリーニング手順が強く望まれる。
❚補足2❚
スピン多重度の異なる電子状態間の非放射遷移。典型的には励起一重項状態から三重項状態への遷移であり、形式的にはスピン禁制である。しかし、スピン軌道相互作用(SOC)は、スピン角運動量と軌道角運動量を混合することで、ISCを可能にする。したがって、ISC効率は、臭素や遷移金属中心などの重原子の存在によって大きく変調される。これらの重原子はSOCを増強し、三重項状態多様体🐾18への状態移動を促進する。
🐾18 三重項状態で構成される空間を指しているようである。原語でも、manifold。
3⃣ 従来の古典及び量子手法における限界
本論文は掲題につき、下記㈠～㈢のようにまとめる。その上で、以下のようにと結論している：｢分子間励起状態モデリングの文脈において、このような限界は、高い精度とスケーラビリティの向上を同時に達成できる、代替的な計算戦略(⇒量子計算🐾19)の必要性を浮き彫りにしている｣。
🐾19 なお、変分量子アルゴリズムについては、｢現実的な光感受性物質モデリングに必要な広い活性空間や完全な吸収特性を扱えない｣と述べている。
㈠ 時間依存密度汎関数法
治療波長域における累積吸収とISC効率の正確なモデリングは、光感受性物質の励起状態の正確な計算に依存する。最も広く用いられている手法である時間依存密度汎関数法(TD-DFT)は、しばしば定量的な精度、場合によっては定性的な精度さえも達成できない。(単一参照法である)TD-DFTは、強い相関または縮退に近い状態を示す系を扱うのが困難である。さらに、TD-DFTの結果は交換相関汎関数の選択に非常に敏感であり、経験的な曖昧さをもたらし、多様な化学系にわたる予測力を制限する。具体例としてTD-DFTは、BODIPYベースの光感受性物質の励起エネルギーに0.3～0.6 eVの系統的誤差をもたらす。これは180～600 nmの波長シフトに相当する。治療波長域の幅が、わずか150 nmであることを考えると、この差は、スクリーニングおよび設計ワークフローの精度を著しく損なう可能性がある。
㈡ 多参照法
EOM-CCSD(運動方程式-CCSD)法や、完全活性空間二次摂動法(CASPT2)などの多参照法は、電子相関をより厳密に捉えることでTD-DFTを改良することができる。多くの系において、これらの手法は励起エネルギー誤差を0.2 eV未満に低減でき、場合によっては中規模系に対して”化学的精度を達成する”ことさえ可能である。
しかし、このような手法の一般的な適用性は、急激な計算スケーリングと大きなメモリフット･プリントによって制約される。より具体的には、EOM-CCSD法は、その基礎となるテンソル縮約のためにO(N6)スケーリングを示し(つまり、比較的小さな分子にしか適用できない)、必要な精度を達成する保証はない。CASPT2法や関連する多参照摂動法では、活性空間は通常約16軌道に制限される(ので、やはり比較的小さな分子にしか適用できない )。
㈢ DMRG法
DMRG(密度行列繰り込み群)法は、従来のCAS(活性空間)ベースの方法に代わる強力な方法であり、はるかに大きな活性空間の取り扱いを可能にする。しかし、基底状態の研究は数百軌道に達しているのに対し、励起状態分光法への応用は、特に広大な活性空間と広範な励起状態多様体の両方が必要な場合、はるかに限定されたままである。多くの一重項状態と三重項状態を同時に記述する必要性、不安定な状態平均最適化、そして動的相関補正のコストが相まって、実用的な分光法研究は、限られた数の活性空間と限られた数の励起状態に限定されてきた。
(2) 光感受性物質の基礎
光感受性物質の治療効果は、2つの重要な特性に依存する。1⃣治療に使用される光に対する感度、2⃣腫瘍組織(がん細胞)を死滅させるROSの生成効率である。
1⃣ 感度
光感受性物質の感度は、治療上重要な近赤外波長域における分子の累積吸収量として定量化できる。PDTの実用例では、この近赤外波長域は約700～850 nmとされている。この波長域の選択は、2つの考慮事項㈠、㈡によって決まる。
使用する光の種類は、㈠本質的に人体に有害であってはならない、㈡深部腫瘍を治療するために人体組織に浸透できるものでなければならない。㈠により、紫外線は除外される。㈡により、近赤外線領域が自然に選択される❚補足3❚。
❚補足3❚
人体組織を構成する2つの主要成分である｢水とヘモグロビン｣の吸収特性が、治療波長域を大きく左右する。水は900 nmを超える波長で自然に吸収が強く、ヘモグロビンは約600～650 nm未満で吸収が強い。このことから、700～850 nmの近赤外線領域が自然に選択される。近赤外線のさらなる利点は、一般に波長が長くなるにつれて光散乱が減少することである。つまり、近赤外線の波長は組織に数ミリメートルから数センチメートル浸透し、可視光よりもかなり深くまで達することができる。したがって、高品質の光増感剤は、この治療波長域で、強い累積吸収を持たなければならない。
2⃣ ROS生成効率
酸素分子は自然に三重項基底状態にある。一方、最も一般的な光感受性物質は、一般に一重項基底状態にある。光吸収だけでは分子のスピン状態は変化しない。光感受性物質が効果を発揮するには、｢励起一重項状態を、励起三重項状態に変換する｣内部プロセスを経る必要がある。
三重項状態になった分子は、近くの分子状酸素にエネルギー若しくは電子を移動させることで、基底状態に戻る。タイプ II 経路では一重項酸素を、タイプ I 経路ではスーパー･オキシド及び関連ラジカルを生成する。この一重項状態から三重項状態への内部遷移は『項間交差(ISC)』として知られており、その効率によって、長寿命三重項状態に到達する分子の数が決定される。多くの光感受性物質ファミリーにとって、ISCの効率が高いことは、ROS 形成を直接促進するため、望ましい特性である。
(3) 為参考：溶媒効果
以下、長々と書いているがポイントは、本論文手法は溶媒効果と相性が良い、ということである。
光感受性及び、効率的なROS生成を表象する❶治療波長域における累積吸収と❷ISC 効率を適切にシミュレートするには、周囲環境の影響を考慮することが不可欠である。血漿や細胞内液などの生物学的環境において、溶媒効果は光感受性物質の光物理的挙動を形作る上で無視できない役割を果たす。これらの効果は励起エネルギーをシフトさせ、遷移双極子モーメントを変化させ、スピン軌道相互作用(SOC)に影響を与える可能性があり、これらはすべて❶及び❷に影響を与える。溶媒誘起効果を無視すると、生物学的条件下での光感受性物質の挙動をモデル化する際に誤差が生じる可能性がある。
PDTでは、水が電子吸収または振動吸収をほとんど示さない700～850 nmの範囲で強い吸収を示す光感受性物質を探索する。これにより系は非共鳴領域となる。このため、周波数に依存しない分極連続体モデル(PCM)を用いて、溶媒効果を溶質ハミルトニアンに対する静的補正として扱うことができる。つまり、溶媒ダイナミクスを明示的にモデリングすることなく、溶媒効果を取り扱うことが正当化される。溶媒誘起シフトは吸収特性の位置に影響を与える可能性があるが、これらの効果は周波数非依存PCMの静的分極を通して捉えられる。このアプローチは本論文の設定に適切であり、溶媒共鳴が電子遷移と直接重なる場合にのみ必要な動的溶媒応答のモデリングに伴う複雑さを回避できる。
(4)BODIPY誘導体
様々な光感受性物質の中でも、ジピロメテン･ホウ素錯体(BODIPY)誘導体は、以下に示すような有利な特徴を多数持つ：高いモル吸光係数、狭い発光波長。生物学的環境において、光化学的及び熱的安定性に優れている。バックグラウンド吸収を最小限に抑えながら、700～850 nmの治療波長域を実現可能。水溶性で生体適合性あり。生物学的環境において、光が存在しない状況下では、最小限の副作用しか示さない。
ハロゲン化、ドナー-アクセプター置換、遷移金属の導入といった戦略は、SOCを増強し、三重項形成を促進することが知られている。ただし、それらの❶治療波長域における累積吸収及び❷ISC効率への影響を予測することは容易ではない。そのため、特定の置換が❶及び❷にどのように影響するかを体系的に評価できるシミュレーションが重要となる。
【3】本論文の技術的要素
(0) 全体概要
0⃣ 振り返りと指針の設定
光感受性物質の治療効果を表象する特性｢光感受性及び、効率的なROS生成｣であった。この2つは、❶治療波長域における累積吸収と ❷ISC 効率、という 2 つの｢量子力学的計算変数｣に置換できる。従来の光感受性物質のほとんどは、❶と❷が劣っていた。つまり、近赤外域ではなく可視光域で自然に吸収が生じる。加えて、スピン軌道相互作用≒ISC効率が低かった。❶は｢光の浸透を制御｣し、吸収を赤方偏移させる。❷は｢反応性三重項多様体への状態移動を制御｣し、三重項収量と ISC効率を高める分子編集を誘導する。結果として、PDTの有効性が高まる。
1⃣ 指針に沿った具体策
❶及び❷を効率的にシミュレーションすることで、❶及び❷を考慮した候補物質をスクリーニングし、より効率的な光増感剤の開発を導くことができる可能性があるのであった。❶と❷を、効率的にシミュレーションするには、光感受性物質における主要な励起･緩和経路をモデル化する必要がある。
具体的に本論文では、3つの誤り耐性量子アルゴリズムを開発した。❶に対する誤り耐性量子アルゴリズムとして、しきい値射影アルゴリズム(☞(1))を開発した。❷に対しては、時間発展プロキシ量子アルゴリズム(☞(2))を開発した。3つ目の誤り耐性量子アルゴリズムは、振電動力学アルゴリズムである(が割愛した)。
(1) しきい値射影アルゴリズム
0⃣ 概要
㈠ 端的に・・・
治療波長域における累積吸収(❶)を評価するアルゴリズムである。基本発想は、❶の評価を量子計算に置換して、当該量子計算を効率的に行う、というもの。具体的には、1⃣量子信号処理(QSP)に基づくスペクトル･フィルターと、2⃣低ランク･テンソルハイパー縮約(THC)因子分解を用いた量子ビット化、を用いる。1⃣にて、効率的な量子計算を実現する。2⃣は、量子回路深度の低減に貢献する。
㈡ 噛み砕いて、技術要素を頭出し
❶は、治療波長域内における基底状態から励起状態への｢電気双極子演算子の(重み付き)遷移確率の総和｣を表す。もう少し平たく言うと、累積吸収は、｢双極子励起状態と、治療波長域へ射影した双極子励起状態との重なり｣の2乗で与えられる。重なりは、二値分類に落とし込まれて、効率的に識別されるような工夫が施されている。操作的に言うと、この重なりは、｢効率的な二値分類を実行する｣射影回路の補助レジスタをサンプリングすることによって推定される。
手続きを述べると、まず、対象量子系の双極子遷移振幅を符号化した初期状態を準備する❚補足4❚。次に、対象となる治療波長域に射影演算子を合成する。本論文では、量子信号処理(QSP)を使って❚補足5❚射影演算子を合成する。QSPを使って、所望のしきい値関数への多項式近似を実現し、より少ない回路深度で鋭いフィルタリングを実現する。
❚補足4❚
初期状態は、十分に高い精度で、古典的に取得できると仮定している。本論文では、この仮定の成立は、DMRG(密度行列繰り込み群)法や選択配置間相互作用(CI)法などの古典手法を使用することで、十分可能とする。なお、電気双極子演算子はユニタリ演算子ではない。このため古典的に用意した初期状態である波動関数は、量子符号化する前に、正規化する必要がある。
状態準備にかかる計算コスト(複雑性)は、O(dlog(d))である。ここで、dは電気双極子演算子(行列)の非ゼロ要素の数である。電気双極子行列は疎行列なので、dは小さい(はずであり、状態準備にかかる計算コストは小さい見込み)。
❚補足5❚
正確には、QSPの制約を取り除いた一般化量子信号処理(GQSP)を使っている。量子アルゴリズム･プリミティブはQSPであることを示すため、｢QSPを使っている｣と表記した。また、量子ビット化(Qubitization)も併用している。量子ビット化を使うことで、浅い回路深度を実現させている。
1⃣ 量子信号処理(QSP)に基づく射影演算子を使ったスペクトル･フィルター
光感受性物質における励起･緩和経路を記述するハミルトニアン(システム･ハミルトニアン)のスペクトル範囲内[最小、最大]で、治療波長域[低、高]＝しきい値、に対応する範囲を分離する｢スペクトル･フィルター｣を構築する。このスペクトル･フィルターは多項式近似されるが、その近似多項式はQSPを使って、効率的に生成される。QSPの実行には、ソフトウェアパッケージpyqspを使用している。なお、単一の狭いしきい値をフィルターとするのではなく、2つのヘビサイド関数でフィルターを実装している(ヘビサイド関数が多項式近似される)。これは、その方が安価であることが分かった、ためである。
システム･ハミルトニアンHの固有状態が、しきい値＝[低、高]内か否かは、 単一量子ビットの読み出しにより、明確に区別可能である。操作上は、固有状態に対応する固有値(固有エネルギー)Eに対して｢Eー低｣と｢Eー高｣の符号を決定することに帰着する。そして、これは期待値⟨E|H′|E⟩の符号を抽出することによって実現される。ここで、H’はHをしきい値でシフトしたハミルトニアンである。
2⃣ 低ランク･テンソルハイパー縮約(THC)因子分解を用いた量子ビット化
本論文では、量子ビット化を用いて、QSP射影演算子の実装を行っている。トロッター公式を用いて実装する方法もあるが、十分に小さな誤差を達成するには、多数のトロッター･ステップが必要となる。つまり、回路深度が深くなる(ので避けた)。
標準的な量子ビット化では、ハミルトニアンはユニタリ行列の線形結合(LCU)形式で表される。1⃣で示したように、❶を評価するアルゴリズムは、期待値⟨E|H′|E⟩の符号を識別することに帰着する。従って、シフトされたハミルトニアンH’のLCU表現が構築される。H’のLCU表現は、低ランクのテンソル･ハイパー縮約(THC)形式で表現される。
(2) 時間発展プロキシ量子アルゴリズム
0⃣ 概要
ISC効率の代理変数を計算することで、ISC効率を効率的に｢推定できる｣アルゴリズムである。ISC効率を正確に予測するには、スピン軌道相互作用(SOC)を高い忠実度で扱う必要がある。しかし、そのために必要な計算量は膨大である。そこで、時間発展プロキシ量子アルゴリズムは、SOCによって誘起される短時間ダイナミクスを解析し、一重項-三重項混合の最低次🐾20の混合を捉える代理変数を計算する。
あくまで代理変数なので、絶対評価はできないが相対評価は可能である。つまり、この代理変数があれば、光感受性物質の候補をISC 効率基準でスクリーニングすることが可能となる。具体的には、既知のISC効率を持つ参照分子に対して代理変数を計算し、参照分子の代理変数値を超える候補物質を拾い上げれば良いことになる。時間発展プロキシ量子アルゴリズムを表す量子回路は、以下の3つに分解できる：1⃣初期状態の準備、2⃣SOCハミルトニアンによる短時間時間発展、3⃣修正アダマール検定による測定。
🐾20 英語で言うとLeading order。最低次という日本語の印象とは裏腹に、最も重要な次数という意味になる。(第0次、)第1次、第2次、･･･と書き下すと、英語と日本語の齟齬が解消する。
1⃣ 初期状態の準備
従来、ISC効率を評価する上での主要な課題の一つは、状態準備にあった。時間発展プロキシ量子アルゴリズムでは、まず、全ての光学活性励起状態の双極子作用重ね合わせを準備する。これらの初期状態には高エネルギー成分が混入するので、低エネルギー部分空間に確率的に射影する演算子を適用して、混入を抑制する🐾21。
🐾21 しきい値射影アルゴリズム(☞【3】(1))を援用する。
2⃣ 時間発展
SOCハミルトニアンHSOCは一体演算子であるため、その下での時間発展はfast-forward可能❚補足6❚である。本論文ではHSOCの対角化に際して、基底変換を適用することで、ブロック符号化の必要性を軽減している。この基底変換は、多体ユニタリ演算子を誘導する。この多体ユニタリ演算子はThoulessの定理🐾22を用いて生成することができ、Givens回転🐾23を用いて分解することができる。
最終的な時間発展演算子は、これらのユニタリ演算子とパウリZ回転の積で表される。
🐾22 真空状態に、指数関数の(肩に乗る)形で｢粒子の生成演算子並びに正孔の生成演算子｣を掛けていくと、真空状態と直交しない全てのスレーター行列式を作ることができる、という定理らしい。尚、そのような｢直交しない全てのスレーター行列式｣を利用すると、スレーター行列式で表される多体波動関数が、多体系の状態全てを含んでいることを保証できる。
🐾23 非対角要素をゼロにできる相似変換。一度に、ゼロにすることができる非対角要素は2つのみであるが、繰り返し適用することで、全ての非対角要素を0 にすることができる。
WD(Δt)NW†＝WD(NΔt)W†＝WD(T)W†＝U(Δt)N≒U(T)となるようなユニタリ演算子の存在は、サウレスの定理によって保証される(という理解で良いのだろう)。対角行列Dを実際に作る場合、ギブンス回転を使う(という理解で良いのだろう)。
3⃣ 測定
標準的なアダマール･テストのセットアップは、{”単独の”入力状態 |ψ⟩、補助量子ビット、ユニタリ演算子U(t) = exp(−itHSOC)に対する制御ゲート}である。アダマール･テストが推定する期待値は、⟨ψ|U(t)|ψ⟩の対角要素となる。
一方、修正アダマール･テストでは、一重項と三重項から選択された2つの異なる参照状態の重ね合わせを、入力状態とする。この入力状態を使用すると、一重項成分と三重項成分の干渉により、非対角行列要素に直接アクセスできるようになる。この修正により、非対角遷移の振幅を直接測定することができる。
【4】推定結果
(0) 全体俯瞰
以下4つ⓵～⓸のBODIPYに対して、❶治療波長域における累積吸収(☞下記(1))、❷ISC効率(☞下記(2))を量子計算するために必要な量子リソースが推定されている：⓵オリジナル🐾24BODIPY、⓶臭素化アザ🐾25-BODIPY、⓷トリアゾリル🐾26置換アザ-BODIPY、⓸白金-BODIPY錯体。
⓵は、以下の理由からベースラインとして扱われる：⓵は、CAS(12,11)つまり電子数12、空間軌道数11の完全活性空間に対して、CASPT2(完全活性空間二次摂動法)計算で、実験とほぼ一致する結果が得られている。従って、古典計算でカバー出来る大きさの分子ということになる。
量子リソースは、論理量子ビット数とトフォリ･ゲート数である。下記(1)及び(2)では、{論理量子ビット数、トフォリ･ゲート数}として表している。Nは、活性空間の空間軌道数である。トフォリ･ゲート数の単位は107であることに注意。
なお、量子リソース推定には、XanaduのPennyLaneが用いられている。
🐾24 オリジナルの意味は、ベースライン(基準)として使うという意味であろう。
🐾25 アザ(aza)とは窒素を意味し、アザ置換とは、炭素環の炭素を窒素で置換することを意味するらしい。
🐾26 トリアゾール環(窒素原子3つを含む五員環)がアザ-BODIPYに結合している、という理解で良いのか、よく分からない。
(1) しきい値射影アルゴリズムの量子リソース推定
⓵→N＝(11,15,19)。N=11:{177、2.72}、N=15:{197、8.44}、N=19:{223、17.9}。
⓶→N＝(17,21,45)。N=17:{211、7.64}、N=21:{228、16.3}、N=45:{346、248.0}。
⓷→N＝(11,19,35)。N=11:{177、2.88}、N=19:{220、12.4}、N=35:{353、97.5}。
⓸→N＝(16,24,30)。N=16:{207、7.10}、N=24:{246、29.1}、N=30:{271、60.5}。
(2) 時間発展プロキシ量子アルゴリズムの量子リソース推定
⓵→N＝(11,15,19)。N=11:{176、7.96}、N=15:{196、24.3}、N=19:{222、51.3}。
⓶→N＝(17,21,45)。N=17:{210、22.0}、N=21:{227、46.7}、N=45:{342、578.0}。
⓷→N＝(11,19,35)。N=11:{176、8.42}、N=19:{219、35.5}、N=35:{299、324.0}。
⓸→N＝(16,24,30)。N=16:{206、20.4}、N=24:{245、83.1}、N=30:{270、173.0}。
【5】考察
(0) PDTの光感受性物質設計に適用される高精度でスケーラブルなアルゴリズムは、量子アルゴリズムになるざるを得ないであろう･･･という論拠は、治療波長域が700～850nmと非常に狭い、という事実に由来する。まずは、そのようなユースケース(アプリケーション)を見つけたことに、賞賛を送るべきであろう。
(1) 本論文手法が、高い精度とスケーラビリティを両立できる訳、改めて整理しよう。まず、PDTの治療効果を、量子力学的計算変数に置換したこと、が上げられる。量子力学的計算変数は２つあり、それぞれに量子アルゴリズムを開発した。しきい値射影アルゴリズムが｢高い精度とスケーラビリティを両立できる｣理由は、しきい値射影アルゴリズムを、量子信号処理(QSP)を量子アルゴリズム･プリミティブとして適用できる枠組みにしたから、と考えられる。ちなみに、QSPで生成する、ヘビサイド関数を近似する多項式の最適な次数は、
321.2051＋4.7571/遷移幅
で与えられている。確かに、300次を超える多項式を(古典的に普通に)生成するのは大変そうである。
時間発展プロキシ量子アルゴリズムが｢高い精度とスケーラビリティを両立できる｣は、以下の合わせ技4本㊀～㊃である。㊀計算コストが高い｢絶対量｣を計算せずに、代理変数として十分な｢相対量｣を計算している。㊁状態準備にかかるコストを、しきい値射影アルゴリズムを使って削減している。㊂fast-forward時間発展を適用している(☞❚補足6❚)。㊃標準アダマール･テストを修正した測定を行っている。
(2) しきい値射影アルゴリズムにおいては、｢論理量子ビット数177、トフォリ･ゲート数 2.72×107｣が、古典計算で対応可能であると解釈できる。時間発展プロキシ量子アルゴリズムにおいては、｢論理量子ビット数176、トフォリ･ゲート数7.96×107｣が、古典計算で対応可能と解釈できる。トフォリ･ゲートをCNOTゲートに置き換える(☞❚補足1❚)と、概ね1.6～4.8×108個のCNOTゲートが、古典計算で対応可能と解釈できる。
つまり、2量子ビット･ゲート数が108程度では、量子優位性(量子有用性)は現れないと推測される。QuOpsでいうと、最低でもギガQuOpsが必要ということになる。確かにそれは、誤り耐性量子コンピューターの領域であろう(し、利用できる時期は、やや先であろう)。
古典コンピューターand/or古典アルゴリズムの進化は止まったわけではないから、ギガQuOpsでも十分とは言えない。テラQuOpsまで行けば、セーフティ･ゾーンであろうか(ただ、本当に実用的なアプリケーション実行に必要とされるのは、テラQuOpsと言われてきたので、違和感はない)。ちなみに、英オックスフォード･クォンタム･サーキッツ(モダリティ：超伝導)は2031年にギガQuOps、34年にテラQuOpsを達成すると予告している(https://oqc.tech/technical-roadmap)。
(1) 背景
米国防総省傘下の国防高等研究計画局(DARPA)が2021年から実施している量子ベンチマーク(QB)というプログラムがある。QBの目的は、｢誤り耐性量子コンピューター(FTQC)の開発競争において、何に重点を置くべきか？を正確に測定できる、量子コンピューティングの指標を開発すること｣である。指標開発に続いて、当該指標を検証可能にし、重要な性能しきい値に到達するために必要な量子リソース及び古典リソースを推定すること、もQBの目的である(☛https://www.darpa.mil/news/2021/utility-quantum-computers)。
QBプログラムのサイト(https://www.darpa.mil/research/programs/quantum-benchmarking)には、11件の研究成果が掲載されている。本稿で扱う研究はその内の一つであり、数値流体力学(CFD)における等温アプリケーションである。
(2) 当初は、量子リソース推定だけだった？・・・
L3Harris🐾1他🐾2の研究者は、非圧縮性CFDにおける等温アプリケーションについて、研究を行った。具体的には、船体設計における抗力計算を扱った。基礎方程式は(量子コンピューターの使用を見据えて)、ナヴィエ･ストークス方程式(NSE)ではなく、格子ボルツマン方程式(LBE)とした。さらに、非線形問題を扱えない量子コンピューターでも取り扱えるように、LBEに対して、カールマン線形化を施す。
当該研究の目的は、｢投資家が、量子処理装置(QPU)開発に資金を提供する価値があるか否か、を判断するために必要な情報を提供すること｣である。アウトプットは量子リソース推定値であり、具体的に言うと、論理量子ビット数とTゲート数を推定した。研究成果は、24年6月(第1版)、24年9月(第2版)そして25年12月8日(第3版)に、arXivにて発表された。本稿は、第3版(以下、本論文[*260])🐾3を扱う。
🐾1 米Harris Corporationと米L3 Technologies社が合併して誕生した企業。2019年7月1日からL3Harris Technologies, Inc.となった。産業セクターは航空宇宙･防衛で、防衛産業部門で売上高米国第6位らしい。
🐾2 米ジョージア･テック･リサーチ･インスティテュート(GTRI)、米ザパタ、米MIT(マサチューセッツ工科大学)リンカーン研究所。GTRIは、米ジョージア工科大学のキャンパス内に設置された米軍の研究機関。MITリンカーン研究所は、米国防総省とMITとが共同で設立した研究所。ザパタは、米ハーバード大学発の量子S/Wスタートアップ。設立は2017年であるが一度、事業停止(24年10月、ザパタAI)。25年9月、ザパタQuantumとしてリブート。
🐾3 ver2までの論文タイトルはFeasibility of accelerating incompressible computational fluid dynamics simulations with fault-tolerant quantum computers(誤り耐性量子コンピュータを用いた非圧縮性数値流体力学シミュレーションにおける計算加速の実現可能性)。ver3は、Detailed assessment of calculating drag force with quantum computers: Explicit time-evolution precludes exponential advantage for nonlinear differential equations(量子コンピューターを使った抗力計算の詳細な評価: 陽解法は、非線形時間発展方程式の指数加速を排除する)。
(3) CFDで指数加速が実現するという研究成果が・・・
arXiv版ver2の第7章結論は、以下のように結ばれている：｢より成熟した量子コンピューティング･アプリケーションで見られる桁違いの改善を考慮すると、FTQCが、将来的に非圧縮性CFDアプリケーションで有用性を発揮するかどうかを見極めるためには、さらなる研究開発が必要である｣。さらなる研究開発の成果は、米パシフィック･ノースウェスト国立研究所他の研究者により、2025年1月10日@Physical Review Research[*170]にて発表された🐾4。これは、CFDで指数加速が実現するという成果である(☛詳しくは、こちらを参照)。本論文は、その成果を否定する。
本論文の技術的工夫は、独自のブロック埋め込み手法を開発したこと及び抗力を推定する量子プロトコルの開発である。しかし本稿では、そちらには焦点を当てず、｢CFDで指数加速が実現できない｣との主張に焦点を当てる。
🐾4 なお、arXivではver1が23年3月、最新版であるver3が24年3月に公開されている。
【1】本論文の主張
本論文は、以下を主張する：
(1) 量子コンピューターを用いたCFDで、指数加速🐾5は期待できない。
※ 実のところ、射程はCFDに限定されず、影響は甚大である。fast-forwardできない時間発展方程式の陽解法では、指数加速が期待できない(☞【4】(2)1及び5⃣⃣)。
(1-1) 一応、エンド･ツー･エンドでの評価である。(☞【3】(1)及び(2))
(1-2) 本論文の射程は陽解法である。しかし、陰解法だからと言って、指数加速は期待できない。(☞【4】(2)4⃣)
🐾5 2次加速に満たない僅かな加速(1.11次加速)は期待できる(☞【4】(3))。ただし、量子誤り訂正に要するオーバーヘッドを鑑みると、それさえ消える可能性がある。
(2) LBEは、量子線形ソルバーを適用する場合でも、法外なコストがかかる。実用的🐾6な問題スケール(Re～108)での計算量は、1039である。
🐾6 本論文はCFDの対象を、船としている。船の場合、代表的なReは108らしい。ただし、たとえ(産業レベルではあり得ないほど小さな)Re～101だとしても、計算量は1021であり、法外であることに変わりはない。
(3) なお本論文には、量子CFDソルバーの市場規模について議論されているが、本稿では割愛した。
【2】事前整理
(0) 数値解析処理のまとめ
本論文と[*170]は、全く同じ処理を行っている。
step1・・・CFD問題の基礎方程式を、NSE(ナヴィエ･ストークス方程式)からLBE(格子ボルツマン方程式)に変換する。
step2・・・衝突項は、一般的なBGK(Bhatnagar-Gross-Krook)形式に従い、3次多項式で近似する。並進する離散速度の方向モデルは、D3Q27🐾7である。
step3・・・カールマン線形化を実行する→カールマン線形化LBEに変換。カールマン線形化によって、有限次元非線形偏微分方程式を、無限次元線形常微分方程式系に変換する。
step4・・・カールマン線形化における打切り次数を、3次とする。
step5・・・(3次)線形常微分方程式系を、線形代数方程式系に変換する。
step6・・・量子線形システム･アルゴリズム(QLSA)を使って、線形代数方程式系を解く。
🐾7 3次元かつ、離散分子速度の数が27であるモデル。
※ [*170]の補足━━━入力問題(☞【3】(1))と出力問題(☞【3】(2))は、無視している。つまり、エンド･ツー･エンドではない。
※ 本論文の補足━━━入力問題と出力問題に対応している。つまり、(一応)エンド･ツー･エンドである。ただし、入力問題を回避するような”テクニカル”な境界条件を課している(☞【3】(1))。このため、流体解析対象の幾何形状は、必ずしも実用的ではないことに注意。出力問題は、真正面から回避している(☞【3】(2))。
(1) 先行研究の発見
[*170]は、カールマン線形化LBEを量子線形ソルバーで解くタスクにおいて、以下を発見した：クエリの複雑さ(2 量子ビット･量子ゲート数とオラクル･クエリの数)が、空間格子ノード数nに対して、poly(log(n))である。古典DNS(直接数値計算)🐾8ではnに対して、ほぼ線形である。これは、量子アプローチが古典アプローチに対して、指数加速🐾9を実現することを示唆している。
🐾8 量子化学計算であればFCI(完全配置間相互作用法)に相当する。
🐾9 n＝exp(log(n))なので、古典:exp(log(n))、量子:poly(log(n))である。poly(多項式)→exp(指数関数)は、指数加速である。
(2) ただし書き
ただし、上記(1)で述べたpoly(log(n))の実現には、次の条件がつく：格子時間ステップTが、空間格子ノード数n に対して、多項式依存性(若しくは、さらに悪い依存性＝指数関数依存性)を持っている場合、指数加速は無効になる。本論文は、T ∝ n1/3である(多項式依存である)ことを示し(☞【4】(2)2⃣及び3⃣)、指数加速が無効になる、と主張する。
【3】量子リソースの推定
(0) 前捌き
0⃣ タスク
船舶を対象とした、抗力を計算するタスクを採用している(抗力計算を選択した理由は、☞下記(2)を参照)。本論文は実用的であることに拘っているので、実在する船舶3隻のデータを使っている。Reは2隻が106のオーダーで、1隻が108のオーダーである。
ただし、入力問題(☞下記(1))を避けるため、形状は実際の船舶のような複雑形状ではなく、単純な幾何形状(球)を仮定し、当該形状を通過する流れの抗力を計算している。
1⃣ 量子リソース推定の対象
具体的な量子リソース推定値は、論理量子ビット数とTゲート数である。論理量子ビット数は、量子計算の実行に必要な空間リソースの代用として機能する。Tゲート操作は、論理ゲート操作における律速操作と考えられている。そのため、実際の量子計算の実行時間は、Tゲート数にほぼ比例すると予想される。故に、T ゲート数は量子計算の実行時間の代用として機能するとして、量子リソース推定値として選択されている。
2⃣ 方法論
量子アルゴリズムがTゲートを呼び出す回数を数える。その回数を、Tゲート数とする。論理量子ビット数は、状態準備、量子振幅推定、ブロック埋め込み、に必要な論理量子ビット数を数えている。
LBEの場合🐾10、輸送(あるいは並進、ストリーミング)行列及び衝突行列の効率的なブロック埋め込みによって、量子リソースが最も大幅に削減される。本論文では、輸送行列及び衝突行列の効率的なブロック埋め込み方法を、独自に開発している(本稿では、扱っていない)。
🐾10 量子化学アプリケーションの場合も、効率的なブロック埋め込みを使用して、量子リソースを大幅に削減できることが示されている。量子化学の場合は、ハミルトニアンの対称性を利用する。
(1) 入力問題の回避：実用的ではない形状
0⃣ 為念：前提知識
量子コンピューターによる計算は煩わしいことに、㊀初期量子状態の準備、㊁計算、㊂測定の3段階から成る。㊀は、古典情報を、量子回路を使って量子状態(量子情報)に符号化する。㊁は、始まりから終わりまで、量子状態のまま行われる。㊂は㊁で得られた計算結果(＝量子状態である)を、(適当な･･･通常は、計算基底を使った)測定によって、古典情報に変換する。
CFDは時間発展問題であり、CFDを解く場合には一般的に(実用的な問題であれば必ず)、初期条件と境界条件を与える必要がある。㊀：古典情報を(量子回路を使って)量子状態に符号化する作業は、多大なコストを要することが知られている。一般的には『入力問題』というピンとこない名称で呼ばれる(あるいは、データ･ローディング問題とも呼ばれる)。古典情報量をnビットとして、入力問題は最悪の場合、exp(n)のコストがかかるとされる。㊀がpoly(n)で可能＋㊁が指数加速を実現＋㊂が多項式時間で測定可能、な場合にのみ、エンド･ツー･エンドの指数加速が実現する。
なお、入力問題の対象は、｢初期条件･境界条件｣の他に、(物理定数などの)パラメータも含まれる。
1⃣ メイン
CFD問題を解くために与えられる初期条件と境界条件は、入力問題を回避するために、ある程度の均一性または幾何学的単純性を持たなければならない。本論文で検討するCFD問題は、均一な初期条件と幾何学的に単純な周期境界条件❚補足1❚を課す。具体的には、初期条件として、一様流体速度ベクトル場を選択する。物体表面の境界条件としてBounce-Back(跳ね返り)条件🐾11を課す。船体の形状は、球面形状並びに、離散的な多面体の集合を選択している。
🐾11 時刻tにおいてノードAとノードBを考え、時刻t+1におけるノードBの速度分布関数f(t+1|B)を与える境界条件。
❚補足1❚
入口/出口境界条件は、考慮されていない。対象システムの解析は、一様な初期速度に対応する均質な局所平衡状態から開始され、3次元で周期的境界境界が課される。この選択は、以下の 3 つの理由による。a)入口/出口で跳ね返る集団に単純に固定速度を課すと、領域内に音圧波アーティファクトが生じる可能性がある。b)より複雑な入口/出口条件は、現在も活発な研究対象となっている。c)入口/出口の速度集団は、シミュレーション前には通常不明である問題のバルク物理と整合する必要がある。
ただし、周期境界条件の選択は、2 つの懸念を引き起こす。1 つ目は、入口/出口条件を外部の強制/駆動項として符号化しないため、量子リソース推定値の数が少なくなることである。2 つ目は、速度場を駆動する入口がないと、シミュレーションはゼロ速度場に向かう傾向があることである。ここで懸念されるのは、最終発展時刻Tにおける速度場が入口駆動速度場よりも低くなり、したがって測定可能な抗力も低くなることである。従って、今後の研究では、高度な非反射入口/出口境界条件を導入する可能性がある。
(2) 出力問題の回避：抗力の推定
0⃣ 前提知識の整理
上記(1)0⃣㊂は、同㊀入力問題と同様に高コストであることが知られており、出力問題と呼ばれる。出力問題は、量子加速を阻害する可能性があることが、広く知られている。
古典CFDは通常、時間発展した流体速度ベクトル場全体及び、その他スカラー場を出力する。その後、様々な後処理解析を通じて、主要な設計性能指標を計算することができる。量子CFDで同様の処理を行う場合、量子トモグラフィーの法外なコストのため、流体速度ベクトル場全体を完全に読み出すと計算上の優位性が失われる。量子CFDで量子優位性を実現するには、少なくとも、関心のあるスカラー量のみを計算し、出力する必要がある。
1⃣ 概要
出力問題を効率的にするために、対象の出力は解ベクトルの低次の多項式関数(この場合は線形)である必要がある。本論文では、解状態の『スカラー関数』のみを抽出することに焦点を当て、具体的には、抗力を計算している。
2⃣ メイン：抗力を推定する量子アルゴリズム
㈠ 概要
本論文は、 抗力を推定する独自の量子アルゴリズムを開発している。障害物に対する抗力は、定常流れに対応する、ある発展時刻における符号化された位相空間密度の線形結合として計算できる。抗力は2つのベクトルの内積(＝スカラー)として表されるため、これら2つのベクトルに比例する振幅を持つ、量子状態を準備できれば、(比例定数が既知であれば)反復量子振幅推定などの、振幅推定のための量子アルゴリズムを用いて、抗力を推定できることになる。
㈡ フローチャート
㊀ 抗力係数を符号化するユニタリ演算子を構築し、｢抗力の線形係数を符号化するベクトルの二乗ノルム｣vを決定する(解析的に可能)。
㊁ 解状態を符号化するユニタリ演算子を構築し、｢解(速度分布関数)の線形係数を符号化するベクトルの二乗ノルム｣fを決定(量子振幅推定を使用)。
㊂ 量子振幅推定アルゴリズムを用いて、二乗振幅aの推定値āを得る。
㊃ v、ｆ、āから、抗力の推定値＝v×f×(2ā－1)を得る。
(3) カールマン線形化が有効であるための制限
カールマン線形化では、非線形性が十分に弱くなければならない。LBEでは、これはマッハ数に比例する。本研究で扱う問題は、マッハ数0.3以下の非圧縮性領域に限定される。本研究はReを101～108の範囲で考慮する。108に相当する流体速度は、約100 m/sである。この速度は水中の音速よりも1桁低いため、すべてのパラメータセットが非圧縮性領域に十分収まる(キャビテーション効果は考慮されない)。
【4】量子加速についての判断
(1) 指数加速が期待できないという判断の根拠
レイノルズ数(Re)を101～108まで振って、カールマン線形化格子ボルツマン方程式(CL-LBE)の量子解法に要する量子リソースを推定した。つまり、｢論理量子ビット数(Ⓛ)と Tゲート数(Ⓣ)｣を推定した。Ⓛ×Ⓣは、計算量の代用として使えるのであった。log(Re)とlog(Ⓛ×Ⓣ)とをグラフ上にプロットしたところ、見事に線形関係が現れた❚補足2❚。その傾きは凡そ、2.68であった。2.68×log(Re)＋β＝log(Ⓛ×Ⓣ) → log(Re2.68)+β＝log(Ⓛ×Ⓣ) → 定数×Re2.68＝Ⓛ×Ⓣ。
まとめると、CL-LBEの量子解法における計算量は、O(Re2.68)ということになる。
❚補足2❚
泥臭いやり方に思えるが、計算複雑性理論に従って計算量のスケーリングを評価する"だけ"なので、このようなやり方で十分ということになる(のだろう)。
(2) 指数加速が期待できない理論的背景
1⃣ メイン
本論文では、指数加速が実現されない理由を、以下のように考えている：時間発展微分方程式系の陽解法によって課される、時間解像度と空間解像度の間の固有の関係。
数値収束のためのクーラン･フリードリヒス･レヴィ(CFL)条件を課すと、時間ステップ数Tは格子点数nに対して、多項式依存性を持つ(☞3⃣)。理論的に予測される指数加速(nに対して、poly(log(n))でスケーリングされる計算量)は、シミュレーション時間(実質的には時間ステップ数T)がnに対して、多項式依存性を持つ場合、打ち消されるのであった。
2⃣ 具体的･定量的に議論する1
本論文では、まず、格子ボルツマン法(LBM)の理論･実践･実装を網羅した入門書[*261]に記されている『非圧縮性領域におけるLBMの安定性と精度を確保するために必要な』関係式を使用している。当該関係式からは、∆t ∝ ∆x2が導かれる。この関係式を基に単純な計算を行うと、多項式依存性T ∝ n1/3が示され、指数加速が無効になる(☞【5】(1))。
3⃣ 具体的･定量的に議論する2
本論文では2⃣とは別に、より一般的な CFL 条件🐾12を使用して議論している。さらに、格子間隔について∆x ∝ 1/Reα 、より具体的には、Reα＝L/Δx(≈n1/3)を仮定している。Lは代表長さである🐾13。この場合も簡単な計算で、多項式依存性T ∝ n1/3が示される🐾14。つまり、2⃣3⃣いずれの場合も、時間ステップ数 Tは全格子点数nのべき乗で増加し、指数加速が無効になる(☞【5】(1))。
🐾12 CFL条件(あるいは、クーラン条件)は、計算の安定性及び、妥当性を保つための必須条件とされる。式で表すと、物理現象における伝播速度×計算時間ステップ≦計算格子幅、である(本論文では、＝が採用されている。つまり、クーラン数＝1である)。
🐾13 船体シミュレーションの場合、Lは船体長である。
🐾14 もう少し丁寧に言うと、T ∝ n(1+α)/3αが示される。n(1+α)/3α＝n1/3α＋1/3なので、αに依らず、多項式依存性(1/3乗)が現れる。なお、1/3乗の『3』は、解析対象の空間次元が『3次元』であることに由来する。
4⃣ 本論文の射程外領域
本論文は掲題に付き、次のようにまとめている：時間ステップ数への線形依存性を回避するには、根本的なアルゴリズムの改良が必要であるか、あるいは、時間ステップ数が空間格子点数に多項式的(あるいはそれ以上)に依存しない応用例に限定する必要がある。例えば、陰解法では一般的にCFL条件の緩和が可能である🐾15。定常状態特性の計算🐾16では、初期条件から時間発展させるのではなく、定常状態を直接解く代替手法を検討する価値があるかもしれない🐾1。また、結合古典振動子のシミュレーションなど、指数加速ではないものの、超2次加速を可能にする微分方程式の応用にも有用である可能性がある。
🐾15 陰解法はクーラン数を大きくとれるだけで、CFL条件を完全に無視できるわけではない。故に、ΔxとΔtとが無関係になるわけではなく、多項式依存性は解消されない。つまり、指数加速は実現しないことになる(はず)。
🐾16 乱流の定常解析には、RANS(Reynolds-Averaged Navier-Stokes Simulation)モデルが用いられる。RANSのまま非定常解析を行うURANSをはじめ、RANSには様々なバリアントが存在する。
5⃣ fast-forwardできないから
また、本論文は次のような解釈を示している：fast-forward(☛こちらを参照)できない流体力学シミュレーションの場合、量子アルゴリズムにおけるクエリ複雑性のべき乗則によるスケーリングが発生する。そのため、指数加速が実現しない🐾17。別の見方をすると、やはり流体解析(乱流解析)と量子コンピューティングの相性が悪い、ということになるだろうか。なお、米Psi Quantum(量子H/Wスタートアップ、モダリティ:光)も、CFDで指数加速が実現しない、と主張している(☛こちらを参照)。
🐾17 fast-forward可能な量子物理シミュレーションや量子化学シミュレーションは、指数加速が期待できる。
(3) 為念：1.11次加速
古典DNSの計算量はO(Re3)とされる。O(Re√3)であれば、2次加速である。前述の通り、CL-LBEの量子解法における計算量は、O(Re2.68)であった。2.68≒30.8973である。1/0.8973≒1.11なので、1.11次加速ということになる。
【5】考察
(1) ⓵【4】(2)2⃣はLBM(LBE)を対象に、文献[*261]から特別な関係式を天下り的に採用した。これだと汎用性が低いので、⓶同3⃣では、時間発展方程式の陽解法に対する基本的な要件である、CFL条件を課した。正確には、∆x ∝ 1/Reαも要請している。これは、レイノルズ数が大きい＝高乱流の場合には、狭い格子間隔を設定するという自然な要請である。⓶の構造を改めて整理すると、{時間ステップ数はΔtと関係あり｜ΔtはCFL条件でΔxと関係あり｜Δxと全格子点数は関係あり}→時間ステップと全格子点数は関係あり＝指数加速なし、となる。
このような自然な要請であるから3⃣の設定は、LBEを越えてNSEに留まらず、時間発展方程式全体を射程に収める、と本論文では考えている。かなり、インパクトは大きい。グーグルの｢問題優先アプローチでは量子優位性の実現は難しい｣という指摘が思い出される。若しくは、諌言耳が痛い･･･
(2) 改めて整理すると本論文が主張する｢指数加速が実現できない｣理由は、計算部分に存在している。つまり悪名高い、入力問題及び出力問題とは無関係である。最終的に線形代数方程式系に変換して解く限り、本論文の主張からは逃れられない、ということになるのだろう。ただし、∆x ∝ 1/Reαが否定できれば、別の展開が開ける。
(3) 本論文における散逸性(非ハミルトン性)への対応は、ブロック埋め込みにより(十分に)対応されているという理解で良いだろうか。
計算流体力学(あるいは数値流体力学：CFD)は、量子コンピューティングのキラー･アプリケーションとして期待されてきたわけだが、25年1月に発表された[*170]は、遂に｢指数加速･実現｣を宣言した。難攻不落の基礎方程式ナヴィエ･ストークス方程式(NSE)では上手く行かないカールマン線形化が、格子ボルツマン方程式(LBE)では魔法のように(?)上手く行くという論文である(☛こちらを参照)。それから約1年後の、25年12月。指数加速は実現しないという論文[*260]が発表された。実は、上手く行かない条件が満たされていたというのである(☛こちらを参照)。そして、同じ25年12月、同じく｢指数関数は実現しない｣ことを主張する論文[*262](以下、本論文)が、Psi Quantum🐾1及びエアバス🐾2の研究者により発表された(25年12月3日@arXiv)。
本論文は、量子コンピューターを使ったCFD(以下、量子CFD)におけるカールマン線形化格子ボルツマン法の適用可能性を、丁寧に検証した論文である(格子ボルツマン法＝格子ボルツマン方程式を使った流体解析)。先行研究の穴(課題)を抽出し、処方箋を提示している。その処方箋を実施した結果として、いくつかの興味深い結果が示されている(ちなみに、論文内で張られている数式のリンクに誤りが存在する)。なお、本論文は基礎編であり、同じ研究者グループによる応用編(☞【7】)が、25年12月5日にarXivで公開されている。基礎編で扱っている流れ問題や境界条件は、｢やや実用性に乏しい｣。【7】応用編で実用性を高めている。
余談ながら、Psi Quantumとエアバスは、26年1月13日、事業提携を発表した(☞【6】(3))。誤り耐性量子コンピューター(FTQC)の航空宇宙分野への応用を推進するらしい。
🐾1 米国の量子H/Wスタートアップ(モダリティ:光)。創業者の豪州人が英国に留学して、米国で起業した。NISQには目もくれず、FTQC構築のみを目指している(マイルストーンは設けず、ゴールのみに向かっている)。
🐾2 言わずと知れた、EUの航空機メーカー。宇宙部門、防衛部門等も併せ持つ。
【1】本論文の主張
以下(0)の制限下で、本論文は(1)～(5)を主張する：
(0) 制限
㊀ 1次元及び2次元流れを対象としているため、3次元流れについて主張を汎化することは困難。
㊁ 扱っているレイノルズ数(Re)が小さい。また、Reの幅(レンジ)が狭いので、主張の汎化には不十分。
㊂ 2次元の場合、カールマン打ち切り次数が低いので、主張の汎化には不十分。
㊃ 流れは、非圧縮性流れに限定。
(1) カールマン線形化手法の性能は、Reだけでなく、対象とする流れの種類に大きく依存する可能性がある。
(2) 流れの種類によって、Reに｢しきい値｣が存在する。Reがしきい値を上回ると、カールマン打ち切り誤差は収束しない(計算が破綻する)。収束挙動は、{流れの種類、Re、カールマン打ち切り次数}に依存する。
(3) 量子CFDが古典CFDよりも効率的になるためには、大きなReが必要である。一方、カールマン打ち切り誤差の収束を保証するために、Re は十分に小さい必要がある。従って、しきい値Reが十分に大きくならない場合、カールマン線形化手法の適用範囲は、制限される可能性がある。
(4) 条件数はカールマン打ち切り次数に依存する。この依存性は、量子CFDが古典CFDよりも効率的になることを妨げる可能性がある。
(5) 量子CFDで、指数加速は実現しない。
【2】事前整理
(0) 先行研究について
【2】では、先行研究における４つの課題(☞(1))と、当該課題に対して本論文が提示する処方箋(☞(2))を整理する。なお、先行研究とは主に、[*170](☛こちらを参照)及び[*260]🐾3(☛こちらを参照)を指している(それ以外の研究は、適宜示す)。
🐾3 正確にいうと、[*260]は一連の作品のver3。本論文の先行研究として言及されている研究は[*260]のver2。
(1) 先行研究における4つの課題
1⃣ カールマン線形化の適用可能性
㈠ 先行研究への批判１・・・不適切な設定による自明な結果
先行研究[*170]では、衝突項のみのD1Q3 モデル🐾4に対して、｢NC = 3 で、誤差🐾5は10−14のオーダーになる｣🐾6とされている。本論文は、この結果に対して次のように批判する：
衝突項のみでは、NC≥ 3 のカールマン打ち切り誤差εCが、正確にゼロであることが解析的に証明できる。NC ≥ 3 の非自明なεCは、衝突項と輸送(ストリーミング)項の相互作用によってのみ生成される。つまり、衝突項のみの数値解析では評価できない。εC ≈ 10−14という値は、カールマン打ち切り誤差ではなく、数値計算スキームの精度によるものと思われる。
🐾4 つまり、1次元モデルである。
🐾5 正確には、相対誤差 = |fCLBE −f|/f、である。fは、格子で離散化された速度分布関数である。上付きのCLBEは、カールマン線形化した格子ボルツマン方程式の結果であることを示している。
🐾6 ちなみに、本論文の表現は、若干異なる。
㈡ 先行研究への批判2・・・おかしな結果
本論文によれば、[*263]は、NC ∈ {1, 2}、Nx🐾7 = Ny ∈ {16,･･･,64}, レイノルズ数Re ∈ [25, 50]というセットアップで、約 10−3という誤差を報告している。そして(これが重要な点であるが･･･)、NC=1の方がNC=2よりも誤差が小さかった。これは、明らかにおかしい。
🐾7 x方向の格子点の数。下付き添え字がyであれば、y方向の格子点の数であり、zについても同様である。
㈢ 結論
本論文は、課題1⃣に対して、次のように結論している：カールマン線形化適用可能性問題は、まだほとんど答えが出ていない。非線形ダイナミクスの線形埋め込みは、あらゆる量子アルゴリズムにとって極めて重要であるため、本論文では、この問題を詳しく検討する。
2⃣ 条件数のスケーリング
㈠ 背景
量子コンピューターを使ってCFDを解くプロトコルにおいて、(どのような手法を採用するにせよ)線形化は欠かせない。量子コンピューターは線形問題しか扱えないからである。線形化された問題に対しては、量子線形ソルバーが適用される。当該プロトコルで(評価すべき)重要な特性は、線形化された問題に現れる係数行列の｢条件数🐾8のスケーリング｣である。
このスケーリングを確認することは、潜在的な量子優位性を主張する上で極めて重要である。なぜなら、不利なスケーリングの下では、量子アルゴリズムは既知の最も優れた古典アルゴリズムよりも大幅に性能が悪くなる可能性があるからである。
🐾8 係数行列に対して定められる｢当該線形問題を、数値的に解くことが適切か否か｣を判定する指標。条件数が小さければ、数値的に解くことが適切となる。
㈡ 先行研究への批判
本論文は、[*170]と[*260]では、条件数のN🐾9依存性が適切に考慮されていないと批判する。具体的な批判は以下の通り：[*170]では、認識誤りに基づく誤った扱いがなされている。[*260]では"仮置き"が行われていて、取り扱いが不十分である。
🐾9 Nは全格子点数。3次元の場合を書き下すと、N＝Nx×Ny×Nzである。Nx等は、🐾7を参照。
㈢ 結論
条件数のスケーリングは、量子加速の議論に重要であり、本論文では、この問題を詳しく検討する。具体的には、NxをReの冪で表すことにより、条件数をReの冪で表現することで、条件数のスケーリングを考えている。
3⃣ 法外に小さな時間ステップ
㈠ 前置き ー 読み飛ばしても良い
DBE(離散ボルツマン方程式)を直接時間離散化する場合、その精度を確保するには 時間ステップ∆tは、∆t ≪ τ が必要である。ここで、物理的な緩和時間τは、大気圧＆室温＆空気の場合、～10−10秒である。このため、典型的な時間スケールである 10−2～10−1秒では、法外に多くの時間ステップ(＝1億～10億ステップ)が必要になる。これが、DBEではなくLBEを選択する理由である。LBEでは、∆t ≫ τ とすることができる。
㈡ 先行研究に対する批判
法外な時間ステップ数を避けるために、㊀｢τを格子単位の緩和時間に置き換える｣場合、㊁｢追加の2次メモリ項を導入する｣必要がある。本論文によると、先行研究[*170]及び[*260]のアプローチは、｢㊀＋㊁のメモリ項を削除する｣ことと同等である。そして残念なことに、このアプローチは、｢2次精度を失い、(実は、回避したはずの)法外な時間ステップ数が要求される｣。
㈢ 結論
先行研究([*170]及び[*260])で提案されている量子アプローチは、実際には競争力がない可能性がある。本論文では、ブロック埋め込みの工夫によって、この時間ステップの問題に対処する。
4⃣ 出力問題
㈠ 概要及び補記
出力問題への対応は以下の通りである：(量子計算の結果＝出力、を含めた)量子状態に符号化されたデータは、潜在的に指数関数的に大きなベクトルである。このため、出力対象は流れ場全体(＝ベクトル)ではなく、特定の特徴量(＝スカラー)を抽出することである。このように本論文では、｢量子情報から古典情報を抽出する｣という表現を使っている。
改めて本論文では、出力問題＝最終時刻Tにおける流体に関する古典情報を、量子状態|y(T)⟩から抽出する問題である。解に相当する量子状態の振幅は、｢カールマン･ベクトルの高次ブロック｣に指数関数的に集中している。このため、情報の抽出にはコストがかかる。そして、このコストは、カールマン打ち切り次数NCと、カールマン･ブロックの次数kに依存する❚補足1❚。なお、簡便のため本稿では、先に言及した｢カールマン･ベクトルの高次ブロック｣を、k次カールマン･ブロックと呼ぶ。
㈡ 先行研究に対する批判1・・・成立しない仮定
[*170]では、|y(T)⟩の高次(k番目の)カールマン･ブロックのスケーリングについて、暗黙的な仮定をおいている。この仮定は、自然には成り立たない。この仮定を成立させることは技術的に可能であるが、そのためには、条件数にN🐾9(再)によるスケーリングを追加する必要があり、量子優位性は消失する。
㈢ 先行研究に対する批判2・・・非線形情報は抽出できない
[*260]は、すべてのカールマン･ブロックから情報を抽出することで、出力問題を回避しようとしている。しかし、このアプローチでは、線形ダイナミクスに関する情報のみが含まれているブロックの情報しか抽出できない。
㈣ 結論
非効率的なデータ抽出問題を軽減できない場合、非線形ダイナミクスのシミュレーションに量子優位性はあり得ない。と結論付ける。本論文では、基礎方程式の再定式化という、かなり大掛かりな方法を導入して、出力問題に対処する。
❚補足1❚
正確に言うと、O(√NNC－k)である。古典の情報抽出コストは、O(N)と考えられるから量子のコストが古典を上回るには、NC－ｋ＝0若しくは1が必要である(1だと2次加速)。2の場合は、√N2＝Nなので、等コストになる。
(2) 課題に対する処方箋
1⃣ カールマン線形化の適用可能性
当該課題への処方箋は、詳細な検討である(☞【4】(1))。詳細な検討を通して、次の2つの疑問に答える：
㊀ カールマン打ち切り誤差εCは収束するか？
㊁ εCが収束すると仮定した場合、実用上適切なカールマン打ち切り次数NCの選択肢は何か？
2⃣ 条件数のスケーリング
当該課題への処方箋も、詳細な検討である(☞【4】(2))。前述の通り、条件数がReの冪を使って表現される。
3⃣ 法外に小さな時間ステップ
時間の離散化を古典側に残し、時間連続 DBEを埋め込む代わりに、2 次精度のLBEを直接埋め込むことで、時間ステップが法外に小さいという問題に対処する。
4⃣ 出力問題
DBEからLBEの導出を再定式化する。具体的には『シフトされた非圧縮性格子ボルツマン方程式』を導入する。この方程式は、非圧縮性流れに限定されるという代償として、u⋆🐾10の大きさ|u⋆|に比例するノルムを持つ状態ベクトルgを用いて、流体を表現できるようにする(ノルムはL2ノルム)。さらに|u⋆|をO(1)とするべく、|u⋆|が十分に小さくなるようにシミュレーション･パラメータを選択する。こうすることで、出力問題に対処する。
🐾10 格子ボルツマン法(LBM)における流速。
【3】本論文の技術的要素
(0) 本論文の貢献
0⃣ 明確な貢献
本論文は、量子CFDにおいて先行研究の課題を抽出･整理して、課題に丁寧に対応している。その結果、(カールマン線形化を適用する)量子CFDの理解が深まり、貴重な知見が蓄積された。それらを、本論文の貢献という括りでまとめた。
1⃣ Reの冪による表現(スケーリング)
パラメータや(CFD的)物理量に対して、Reの冪による表現式を提示している。流体屋目線で言うと、非常に便利であるし、使い勝手も抜群に良くなっていると思われる。エアバスの知見であろうか。Reのしきい値発見(☞2⃣)も、このスケーリングのお陰である。以下、表現式を具体的にあげる。ちなみに、β≐3/4、⎾･⏋はガウス記号(天井関数)である：
① Nx＝⎾Reβ⏋
② T＝⎾Reβ(D/2＋1)⏋
③ τ＝1/2＋3/Reβ(D/2－1)＋1
④ u＝1/ReβD/2
①＝🐾7参照。②＝格子ボルツマン法(LBM)におけるタイムステップ総数。③LBM法における緩和時間。④＝🐾10参照。
2⃣ しきい値レイノルズ数ReTの発見
㈠ 総論
数値シミュレーションの結果、しきい値レイノルズ数ReTを発見した。レイノルズ数がReTより大きい場合、カールマン線形化は収束しない。この｢収束しない｣の意味は、カールマン打ち切り次数NCの増加が、カールマン打ち切り誤差εCを減少させるのではなく、むしろ増加させる、という意味である。
㈡ 各論1：1次元
1次元の場合ReT ～102である。
㈢ 各論2：2次元
2次元の場合は流れの種類によって、ReTは大きく異なる結果が示された。テイラー･グリーン渦❚補足2❚の場合、ReTは観察されなかった。一方、ガウス渦双極子❚補足3❚の場合は、ReT ≈ 400と推測された。本論文では、(この差異は)｢カールマン線形化法の性能が、レイノルズ数だけでなく、流れの種類に大きく依存する可能性があることを示唆している｣と解釈している。
㈣ 各論3：3次元
本論文では3次元の数値シミュレーションは行われていないので、3次元の場合のReTは、｢2次元の場合よりも大きいであろう｣ことが、”定性的に”推測されているに過ぎない。
❚補足2❚
周期境界条件下におけるナヴィエ･ストークス方程式の特殊解。初期条件(2次元の場合)は、三角関数で(u(0)＝sin(x)･cos(y)、v(0)＝－cos(x)･sin(y))のようにシンプルに与えられるが、時間発展は複雑である。大域解の存在は不明のはずであるが、解析解が存在する。隣り合う渦の向きは反対であり、粘性散逸を釣り合わせるような力がないため、時間とともに指数関数的に減衰し消滅する。
❚補足3❚
正確な訳語はないと思われる(その時点でポピュラーではない?)が、本稿では｢ガウス渦双極子｣とした。英語では、Gaussian vortex dipole。向きが逆で、強度が等しい2つの渦対構造。ガウス関数型の流れ関数が与えられる。解析解は存在しない。こちらも、(粘性散逸を釣り合わせるような力がないため、)時間とともに指数関数的に減衰し消滅するが、渦の逆回転により自己移流が生じる。このため、テイラー･グリーン渦に比べて、より非線形である。流体力学というよりも、光学で取り上げられることが多いようである。流体力学(CFD)で、取り上げる価値が意味があるか不明。
3⃣ ゲート･コスト
ゲート･コスト(Tゲート数)を、カールマン打ち切り次数とゲート誤差の関数として表現している(☞【5】(2))。敢えて、Reの関数とはしていない。
(2) 独自の技術要素
量子線形ソルバーは、[*264]で提示されているアルゴリズムを使っているので、オリジナルではない。
1⃣ シフト非圧縮性LBE
LBE状態ベクトルをシフトさせた新しいベクトル値関数gを解とするLBEを、本論文ではシフトLBEと呼んでいる。シフトとは、ゼロ速度平衡状態からズレているという意味である。また本論文の対象🐾11は、非圧縮性流体に限定されているので、シフト非圧縮性LBEと呼ばれる。シフトされた変数gを導入することで、再スケーリング手法のように条件数に直接的な悪影響を与えることなくgのノルムを制御でき、出力問題における主要な障害を克服することができる(らしい)。
🐾11 [*170]も[*260]も、非圧縮性(正確には、弱圧縮性)を仮定している。
2⃣ ブロック埋め込み
非ユニタリ性に対処する方法は主に、3つある。最もポピュラーな手法は、ユニタリ拡張(Unitary Dilation)あるいはブロック埋め込みと呼ばれる方法である。非ユニタリ系をユニタリ系の部分系と考えて、非ユニタリ系をユニタリ系に埋め込むという手法である。本論文において、どのようなブロック埋め込み手法が開発されたかは、割愛する([*260]を扱った別稿も、ブロック埋め込みは割愛した)。
3⃣ 出力問題ー抗力推定
本論文は(も)、一応エンドツーエンドを標榜しているので、出力問題に正面から向き合っている。入力問題は、ややスカシ気味であるが【7】応用編で、少し踏み込んでいる(ただし、古典情報を量子情報に効率的に符号化する"純粋な"入力問題に踏み込んでいるわけではない)。出力問題に正面から向き合っているとは言え、基本処理は[*260]と同じである。出力も同じ抗力である(故に、独自性はない)。ただし、LBEベクトルではなく、シフトLBEベクトル(☞【3】(2)1⃣)を用いているため、若干の変更が加えられている。
【4】数値シミュレーションによる評価結果
数値シミュレーションは、1次元(D1Q3モデル)と2次元(D2Q9モデル)に対して実施された。3次元に対しては、計算コストが高いとして、実施されていない。
(1) カールマン打ち切り誤差解析
0⃣ カールマン打ち切り誤差の評価指標
本論文は、カールマン打切り誤差を、RMSE(自乗平均平方根誤差)を使って評価しない🐾12,❚補足4❚。代わりに、次の値を用いる。
評価指標＝max∑[㊀/㊁]/N
㊀＝‖カールマン線形化LBEの解ーオリジナルLBEの解‖
㊁＝対象とする流れ問題の代表速度
㊀が、生の｢カールマン打ち切り誤差｣ということになる🐾13。maxは、全ての時間ステップに対して見つける。∑(総和)は、全格子点に渡って行う。Nは、全格子点の数である。代表速度を導入したのは、Reを使う伏線である。ここでは、カールマン打ち切り誤差をReの冪で表現したいので、代表速度＝ReβD/2とする。ここで、β＝3/4🐾14、D＝次元数(つまり、1あるいは2)である。最終的に、 評価指標は
評価指標＝Re－βD/2max∑[㊀]/N
となる。
🐾12 その理由は、｢異なるRe に対して、離散化パラメータ ∆x と ∆t の値が変化するため｣と説明されている。
🐾13 ノルムは、L2ノルム。
🐾14 古典CFDはDNSを想定しているから、量子CFDはコルモゴロフ･スケールに対応する必要がある。コルモゴロフ･スケールの解像を保証にするには、Re(3/4)×Dが必要とされる。つまり、本論文でβ＝3/4を設定することは合理的である。
❚補足4❚
以前の研究との比較を容易にするために、付録Eで、RMSEを使った結果も示されている(ちなみに、紙面の都合上、図は付録F以降の頁に掲載されている)。1次元の場合は定性的には同じ。2次元の場合も、テイラー･グリーン渦の場合は、定性的に同じ。しかし、ガウス渦双極子の場合は、定性的にも異なる。
1⃣ 1次元の場合における収束
初期条件として、2つパターンを考えている。一つは、正弦波を与えるパターンである。もう一つは、議論をシンプルにするため割愛した🐾15。収束は、ReとNCの2変数で評価される。以下、(NC｜Re)という形式で、使用された値を示す。Re＝30～1000🐾16で、NCは1～5である。(1｜1000)、(２｜200)、(2｜1000)、(3｜100)、(3｜500)、(4｜30)、(4｜150)、(5｜50)
【3】(0)1⃣でも示したように、Reにはしきい値が観察される。収束するか否かは、Reとしきい値の関係によって変わる。Reがしきい値より小さい場合は、カールマン打ち切り次数NCに関わらず、収束する。ただし、NCが大きいほど、速く収束し、(カールマン打ち切り)誤差も小さくなる🐾17。
収束すれば良いという目線であれば、1次元の場合、NC＝2で十分であろう。
🐾15 x方向の場所(位置)に応じて、3つの値を与えるパターンである。ただし、3つとは｢適当な初期値、0、－1×(適当な初期値)｣である。
🐾16 乱流の目安がRe＝2000なので、(正面切っては)乱流を扱っていないことになる。
🐾17 評価指標の数値の大きさ自体に、意味はないだろうから、具体的な数値は示さない。
2⃣ 2次元の場合
初期値は、テイラー･グリーン渦の場合、{u(0)＝sin(x)･cos(y)、v(0)＝－cos(x)･sin(y)}の形式で与える🐾18。ガウス渦双極子の場合は、{u(0)＝∂ψ/∂y、v(0)＝－∂ψ/∂x}の形式で与える🐾19。
1次元の場合と同じように、ReとNCについて、(NC｜Re)という形式で、使用された値を示す。Re＝20～250🐾20で、NC＝1～3である。(1｜100)、(1｜250)、(2｜20)、(2｜150)、(3｜20)
🐾18 実際は、xとyの部分が、x→－π＋2πx/Nx、y→－π＋2πx/Nyとなる。x、yは、それぞれの格子上の位置である。Nxは、ｘ方向の格子点数。Nyは、y方向の格子点数。
🐾19 流れ関数は、exp[－((x－s)2＋y2)/2σ]2ーexp[－((x＋s)2＋y2)/2σ2]。sは渦間の距離。σは渦の幅。
🐾20 おもちゃレベルのレイノルズ数である。
㈠ テイラー･グリーン渦
【3】(0)1⃣でも示したように、この場合、しきい値は観察されない。Reの増加に伴い、誤差が逓減していく。NCが大きくなるに従い、誤差は小さくなる。挙動には変化は見られない。つまり、NCが大きくと速く収束するような挙動は見られない。ただし、NC＝3の場合は、Re＝20のみのシミュレーションであり、データが少なすぎるように思われる。
㈡ ガウス渦双極子
こちらも【3】(0)1⃣でも示したように、テイラー･グリーン渦とは異なり、しきい値が観察できる。特徴的な挙動は、｢Reが増加すると、NC＝2の場合の誤差｣は、｢NC＝1の場合の誤差｣、に急速に近づくことである(☞4⃣)。
3⃣ 3次元の場合における推測
しきい値が、(2次元の場合よりも)大きくなる可能性がある、と推測している。2次元の場合、400程度と推測される(☞【3】(0)2⃣㈢)から、1000程度になる可能性もあるだろう。
4⃣ 重要な指摘
本論文は、次のように指摘する：テイラー･グリーン渦とガウス渦双極子における誤差挙動の違いは、両者の根本的に異なるダイナミクスに起因している可能性がある(☞❚補足2❚及び❚補足3❚)。(そして)これは、『カールマン線形化法の性能が、レイノルズ数だけでなく、対象とする流れの種類に大きく依存する可能性がある』ことを示唆している。
(2) 条件数解析
0⃣ 前説
有限次元非線形偏微分方程式(LBE)→無限次元線形常微分方程式系(カールマン線形化)→有限次元線形常微分方程式系(カールマン次数の打ち切り)→線形代数方程式系(離散化)：係数行列A×未知である解ベクトルx＝既知のベクトルb。条件数κ(A)＝Aの最大固有値/Aの最小固有値、である🐾21。
ここでも、条件数をReを使って表現している。Reを振ってシミュレーションを行う条件数を計算する。Re(の対数)と条件数(の対数)をプロットすると線形関係が確認できた(スケーリングの確認は、この泥臭い作業で十分である)。すなわち、
ln κ(A)= χ･ln Re＋ln c
が確認された。変形すると、以下のようにスケーリングされる。
κ(A)＝ c×Reχ
パラメータ χと定数cは、問題の次元D とカールマン打ち切り次数 NCの両方に依存する。
🐾21 余談：ランチョス法を使って、固有値を求めている。
1⃣ 結果と考察
本論文で示されている数値はχ及びcについて、NCが1次元の場合{1～4}、2次元の場合{1、2}である。次元及びNCの増加と共に、χが増加している。
次元とともに χ が増加することは、古典との比較において、大きな問題とはならない。なぜなら、古典計算量もD とともに増加するからである。しかし、NCへの依存は問題となる。小さい誤差で出力を得るためには、高いNCが必要になる可能性がある。この場合、高いNC→大きなχ→大きな条件数となり、悪性の行列(悪条件)となる。そのため、直接的に計算コストに反映される可能性がある。すなわち、古典との比較において、計算高速化の消失をもたらす可能性がある(☞【5】(1)3⃣)。
【5】推定結果
(1) 量子加速の推定
1⃣ 3次元の場合
｢‖g‖＝O(1)🐾22を要請し、Nx＝Reβを仮定し、チャップマン･エンスコグの式を利用する｣と、量子CFDの計算複雑性は、O(Reβ(D/2＋１))であることが導かれる。ここでgは、シフトLBEベクトルである(☞【3】(2)1⃣)。なお、この場合の計算複雑性は、量子クエリ複雑性であり｢Tゲートの数｣ではないことに注意🐾23。古典CFD(直接数値シミュレーション)の計算複雑性は、O(Reβ(D＋１))であることが知られている。
コルモゴロフ･スケールを解像するβ＝3/4を使うと、3次元(D＝3)の場合で、量子CFDはO(Re15/8)＝O(Re1.875)となる。およそ1.75次加速となり、指数加速どころか、2次加速にもやや満たない。
実は、Nx＝Lx/Δxであり、これはLx/Δx＝Reβを仮定していることに相当する。[*260]で示されたように、そのように仮定すると｢格子点数と時間ステップ間に多項式依存性が発生｣し、指数加速が消える。整合性が取れていると言えるだろう。
🐾22 解状態を得るためのコストは、時間離散化後に得られる線形方程式系の条件数に大きく依存するが、ノルム特性にも依存する。‖g‖＝O(1)は、解状態のノルムが有界であることを表している。
🐾23 従って、Tゲートの数を基に量子CFDの計算複雑性を算出している[*260]との比較はできない。[*260]の結果は、O(Re2.68)であった。
2⃣ 1次元と2次元の場合
1⃣と同じことを1次元と2次元に対して行うと、1次元ではO(Re1.125)、2次元ではO(Re1.5)となる🐾24。
🐾24 D=1では1.5(＝3/4×(D+1))ー0.375(＝3/4×D/2)＝1.125。D＝2では、2.25(＝3/4×(D+1))ー0.75(＝3/4×D/2)＝1.5。
3⃣ 1次元と2次元の場合▪別バージョン
本論文では、別のやり方でも、量子CFDの計算複雑性を算出している。具体的には、｢条件数スケーリング(Reの冪で表現した条件数の計算複雑性)O(Reχ)｣を使って、量子CFDの計算複雑性が計算されている。この計算複雑性は上限値を与えるという建付けなので、｢最悪でも、達成される量子加速｣を表す。
D=1の場合は、NC＝{1,･･･,4}に対して、{O(Re1.167)、O(Re1.691)、O(Re2.283)、O(Re2.792)}となる。D＝2の場合は、NC＝{1,2}に対して、
{O(Re1.588)、O(Re1.936)}となる。
4⃣ 量子加速の判別
本論文では、㊀｢3⃣の量子CFDの計算複雑性｣と、㊁｢古典CFDの計算複雑性｣とを比較することで、量子加速の有無について、判別が行われている(ただし、この判別は量子演算のみで、エンドツーエンドではない。エンドツーエンドは6⃣を参照)。㊀ー㊁＞0なら量子加速はない。
D=1の場合は、NC＝{1,･･･,4}に対してλ＝{0.333、－0.191、－0.783、－1.292}🐾25となる。つまり、NC=1の場合にのみ量子加速が発生し得る。D＝2の場合は、NC＝{1,2}に対して、λ＝{0.662、0.314}🐾26となる。つまり、NC=1及び2の場合に、量子加速が発生し得る。
🐾25 D=1における古典CFDの計算複雑性1.5から、(0.375の代わりに)χの値{1.167、1.691、2.283、2.792}を引くと算出される。
🐾26 D=2における古典CFDの計算複雑性2.25から、(0.75の代わりに)χの値{1.588、1.936}を引くと算出される。
5⃣ 外挿を行って、4⃣を拡張すると･･･
条件数スケーリングかつD＝1の場合、NC＝1しか量子高速化が得られない(負号は、減速を意味する)。D＝2の場合、NC＝2までしか本論文に記載ないので、NC＝5まで外挿❚補足5❚した。その結果は、NC＝2までしか量加速は発生しない、というものであった。
❚補足5❚
条件数の対数とReの対数には、線形関係が成立している(☞【4】(2)0⃣)ので、単純な外挿を実施した。つまり、NC＝2とNC＝1の差分を、順に足し上げてNC＝5までの値を外挿した。
6⃣ 出力推定の計算複雑性
最終時刻における境界壁上の抗力を推定するための計算コスト(計算複雑性)について、結果だけを簡単にまとめる。当該コストは、O(√Nx)となる(Nxは🐾7を参照)。Nxの平方根が現れる理由は、｢出力である抗力の推定を[*260]に倣って、振幅推定のための量子アルゴリズムを用いているから｣である(と理解している)。Nxは、Reを使って、Nx＝L/Δx＝Reβとスケーリングする(と、CFDのドメイン知識を使って合理的に仮定される)。つまり、当該コストは、O(Reβ/2)となる。コルモゴロフ･スケールの解像を可能とするβ＝3/4を設定すると、最終的に、O(Re3/8)が得られる。
2次元、カールマン打ち切り次数＝2とした場合に、条件数スケーリングから示唆される、量子CFDの計算複雑性はO(Re1.936)であった(☞3⃣)。O(Re3/8)を足すと、O(Re2.311)になる。これは古典CFD(DNS：直接数値計算)の計算複雑性O(Re2.25)よりも大きい。故に、量子減速となる。
(2) ゲート･コスト
Tゲート数で計測する。このコストは明らかに、Reあるいは、最終時刻における解状態を取得する成功確率を制御するパラメータWに依存する。依存するけれどもReとWで表現せずに、本論文では、カールマン打ち切り次数NCとゲート近似誤差εの関数として表現している。定性的に言うと、 NCに対して3 次依存性を持ち、１/εに対して対数依存性を持つ。1次元と2次元に対して、定量的な結果を示すと、以下のようになる。εは10－6とした。
1⃣ 1次元
(NC＝1、2、3、4、5)→(19万、69万、167万、329万、571万)
2⃣ 2次元
(同)→(2,260万、8,308万、1億9,922万、3億8,957万、6億7,315万)
3⃣ 先行研究との比較
[*260]では、全格子点数Nと誤差の関数として表現されている。本論文と平仄を併せて(D＝２、N＝102、NC＝3、ε＝10－6)計算すると、1,877万となる(はず)。1億9,922万と比べて、1桁少ない(ただし、比較対象が全く同じというわけではないかもしれない)。
【6】考察
(1) 本論文の目的は、量子CFDにおける｢指数加速は、ない｣ことを所与として、わずかでも量子加速が存在するか？を、可能な限り丁寧に検証することである(という理解)。エンド･ツー･エンドを標榜しているが、入力問題(正確には入力問題に付随する問題)は、やや消化不良である(境界条件に構造(対称性)を付与している)。本論文の著者もそれは認識しており、別論文[*266]を書いている。結論として、入力問題をよりシビアにしても(より現実に近づけても)、わずかな量子加速は実現する、と主張している(ただしReが小さすぎて、まともな議論にならないと思われる。結論として、同意できない)。
アーロンソン･アンバイニスの予想から、データに構造がない場合、指数加速は期待できない。量子CFDの場合、構造の有無に関わらず、一般的には指数加速は実現しない(と考えられる)。
2次元かつ最悪の場合、エンドツーエンドでO(Re2.311)であった(☞【5】6⃣)。そして、これは古典CFDより遅いのであった。2次元だとNC=1で、O(Re1.963)となり、かろうじて、量子加速が発生する。これは、古典CFDで24時間かかる計算が、241.963/2.25＝約16時間で終わることになる。現場目線では、これでも十分な気もするが･･･☞(4)。
(2) 【5】(1)5⃣では、2次元のケースで(カールマン打ち切り次数に対して)外挿を行った。ここでは、
3次元の外挿を行おう。本論文では、3次元の計算は全く行われていない。無理やりNC＝5まで外挿して❚補足6❚計算すると、NC＝5まで、量子加速が発生し得るという結果となった。NC＝5でO(Re2.925)。
ただし、上記は量子演算のみの結果。出力推定まで含めると、NC＝3までしか量子加速は発生しない。NC＝3の場合でO(Re2.755)。これは、古典CFDでの24時間計算を約18時間30分に短縮することに相当する。
❚補足6❚
以下、やや煩雑に見えるだろうが、基本発想は｢D=1とD=2との関係｣が｢D=2とD=3との関係｣にも成立するとの仮定である。具体的には、D=3の場合は、2系列の平均average{㊀、㊁}とした。㊀は、同じNCのD＝2の値とD＝1の値の比をとり、その比をD＝2の値にかけてD=3の値とした(つまりDに関して、等比を仮定した)。㊁では、まず、この比に対して、NCにわたり平均をとった。次に、この平均比をD＝2の値にかけてD=3の値とした。
(3) 前述(☞【0】)の通り、Psi Quantumとエアバスは26年1月13日、事業提携を発表した[*265]。一見すると、量子CFDの未来は明るくないように見えるが、まだ探索すべき箇所は残っている、ということであろうか。何と言っても、量子加速の余地が広い、3次元の解析は未着手である。さらに、｢流れの種類によっては、カールマン線形化の性能が異なるかもしれない｣という気付きが、希望になっているのかもしれない。
また、カールマン線形化に拘る必要もないだろう。クープマン･フォン･ノイマン線形化を適用しても指数加速は無理だろうが、スケーリングが改善する可能性はあるだろう。
(4) (CFDでは有名な)東北大学流体科学研究所と豪メルボルン大学の研究者は、AIでURANS(Unsteady Reynolds-Averaged Navier-Stokes：非定常レイノルズ平均ナヴィエ･ストークス方程式)モデルを構築した、と発表した(26年1月29日)[*265]。DNSに対して誤差が1.27%、計算コストが0.55%とする。24時間×0.55%≒24時間/180＝24×60(分/時)/180≒8分。
もちろん、DNSとRANSを同列で比較することには、無理がある。ただ、[*260]は量子CFDの活躍の場として、定常解析をあげていた。その場所も、失うことになりそうである。
1⃣ セットアップ
同じ研究グループは、[*262]で開発したカールマン線形化LBM(CL-LBM)に対して"拡張"を必要とする3つの異なる流れに、当該CL-LBMが適用可能か?を評価した論文[*266]も発表している(25年12月5日@arXiv)。
3つの異なる流れとは、(1)外力で強制駆動されたテイラー･グリーン渦、(2)蓋駆動キャビティ流れ、(3)円柱周りの流れ、である。(2)は面妖に聞こえるかもしれないが、(3)と共に、代表的な｢ベンチマーク流れ｣である。この(1)～(3)は、｢CL-LBMを、産業上重要な実世界の問題に結び付けるには、⓵非自明な境界条件、と⓶外力による強制駆動、を組み込む必要がある｣という思想に基づいて選択されている。
2⃣ 判断軸
適用可能性の判断軸は、次の3つである：❶"拡張"に伴って計算コストが、どの程度増加するか、❷カールマン打切り誤差εCの誤差挙動、❸古典LBMとのL2相対誤差εL2の誤差挙動。
3⃣ ⓵非自明な境界条件と、⓶外力による強制駆動、に対応する"拡張"の中身
⓵➡ [*262]では、周期境界条件(横文字を使うと、ラップアラウンド条件)という玩具の条件を設定していた。現実的なCFDにおいては、｢壁、流れの入口、流れの出口｣に対応できる、滑りなし条件(→ディリクレ条件)、自由すべり条件(→ノイマン条件)が設定できなければ、話にならない。
⓶➡ 実用的な流れ問題では通常、外力(body-force：体積力)が作用しているから、外力を考慮する。具体的には、(1)を考える。
拡張の中身をシンプルに述べると、ブロック埋め込み方法及び、状態準備ユニタリ量子回路が変わる。ブロック埋め込みが変わり、線形代数方程式系の係数行列が変わるため、条件数も変わる。条件数の変化は、計算複雑性に影響する可能性がある(☞【4】(2)及び【5】(1)3⃣)。状態準備回路が変われば、状態準備コストも変わる可能性がある。
(1) 外力で強制駆動されたテイラー･グリーン渦
既述通り、⓶を考慮するケースである。Re＝50である。実用性に拘っているが、Reは非現実的に小さい(☞(5)補記)。はっきり言って実用性は皆無である。
❶ シンプルに、Reの冪として書き下せば、1.836が推定されている(☞(4)考察1⃣)。
❷ εCが、NC＝1,2,3について検討されている。Reの増加と共に減少する。NCの昇順(＝1→2→3)で、εCの水準(εCを結ぶ曲線の位置)が低減する。NC=3の場合、εCの水準自体が大幅に低く、(曲線の)勾配も急峻である。
❸ εL2が、NC＝1,2について検討されている。NC=1のεL2は、NC=2のεL2より小さい(つまり、結果は奇妙である)。ただ、両者とも、一定の有限値に収束している。
(2) 蓋駆動キャビティ流れ
⓵を考慮するケースである。具体的には、滑りなし条件(→ディリクレ条件)を適用する。Re＝25である･･･(☞(5)補記)。実用性は皆無である。
❶ シンプルに、Reの冪として書き下せば、1.712が推定されている(☞(4)考察1⃣)。
❷ εCは、NC＝1,2について検討されている。Reの増加と共に『増加』している🐾27(☞(4)考察2⃣)。NCの昇順(＝1→2)で、εCの水準(εCを結ぶ曲線の位置)が低減する。
🐾27 この理由について、次のように推測している：蓋の近くで強い速度勾配が発生するためだと考えられる。
❸ εL2が、NC＝1,2について検討されている。NC=1のεL2は、NC=2のεL2より大きい(つまり、結果は真っ当である)。両者とも、一定の有限値に収束している。
(3) 円柱周りの流れ
⓵を考慮するケースである。具体的には、自由滑り条件(→ノイマン条件)を適用する🐾28。Re ≤ 10である･･･(☞(5)補記)。実用性は皆無である。
🐾28 正確には、垂直方向＝周期的境界条件。入口＝滑りなし境界条件。出口＝自由滑り条件(→ノイマン条件)である。
❶ 推定されていない。
❷ εCは、NC＝1,2について検討されている。見たままをありにままに言うと、NC＝1はReの増加で減少している。NC＝2は増加したあと、一定値に収束している。しかし、Reの範囲が小さすぎて、確定的なことは言えない、というのが正しい解釈だろう。なお、NCの昇順(＝1→2)で、εCの水準(εCを結ぶ曲線の位置)が低減する。
❸ εL2が、NC＝1,2について検討されている。NC=1のεL2は、NC=2のεL2より大きい(つまり、結果は真っ当である)。そして両者とも、εL2が増大していく。NC=1の場合は『L2相対誤差が約15%』🐾29で飽和する。NC=2の場合は、｢収束挙動を確認するのに十分な時間シミュレーションを実行できなかった｣と説明されている(☞(4)3⃣)。
🐾29 この誤差水準だと解析は、ヲワっている。
(4) 考察
1⃣ ⓵及び⓶を考慮しない場合におけるReの冪は、1.784と推定されている。(1)の冪1.836及び(2)の冪1.712は、この数値と大きく変わらない、とまず主張する。そして、抗力抽出の計算コストO(Re3/8)を控除しても、2次元古典CFDの計算コストO(Re2.25)よりも小さい(わずかな量子加速あり)と主張している🐾30。
疑問1：⓵及び⓶を考慮しない場合より、(2)の方が計算コストが低いのはなぜか？
疑問2：[*262]ではD=2、NC＝2、Re∈[20,150]についてO(Re1.936)が推定されている。これはReの水準で変わっている、という説明が可能だろう。どちらにしても[*266]のReの水準は小さすぎて、その結果で実用性を議論することは難しいだろう。
🐾30 (1)の場合でO(Re2.211)、(2)の場合でO(Re2.087)。これは、24時間が、(1)約22時間43分、(2)約19時間に短縮することに相当する。
2⃣ 物理的に整合しているため正しいシミュレーションが行われた、と判断されている。そして❸を併せて、｢滑りなし境界条件｣への拡張が可能、と判断されている。
3⃣ 自由滑り条件への拡張可能性判断は、ペンディングということになる。
(5) 補記
(1)のような流れ問題Ⓧ2次元では、Re＝50を超えると、カールマン線形化LBMは破綻するらしい(!!)。それを回避する技術はあるが、計算速度が落ちる[*267]。いずれにしても、現状では量子CFDに勝ち目はない。
【尾註】
*241 Daniele Loco et al.、Practical protein-pocket hydration-site prediction for drug discovery on a quantum computer、https://arxiv.org/pdf/2512.08390
*242 Mauro D'Arcangelo et al.、Leveraging analog quantum computing with neutral atoms for solvent configuration prediction in drug discovery、https://journals.aps.org/prresearch/pdf/10.1103/PhysRevResearch.6.043020
*243 吉田紀生他、3D-RISM/RISM理論の概略 、https://www.jstage.jst.go.jp/article/biophys/51/5/51_5_222/_supplement/_download/51_222_1.pdf
*244 大場正春 、解説｜積分方程式理論による溶液熱力学量の計算 、Netsu Sokutei 30（2）66-75 、2003、https://www.netsu.org/JSCTANetsuSokutei/pdfs/30/30-2-66.pdf
*245 水野章敏 、ホッとひといき｜高温液体の中の秩序 、応用物理 第86巻 第2号(2017)、pp.148-150、 https://www.jstage.jst.go.jp/article/oubutsu/86/2/86_148/_pdf/-char/ja
*246 佐藤航地、(博士論文)深層学習を用いた蛋白質の水和構造予測、https://koara.lib.keio.ac.jp/xoonips/modules/xoonips/download.php/KO50002002-20246413-0002.pdf?file_id=187966
*247 杉田昌岳他、創薬を見据えた3D-RISM理論の応用研究 、生物物理57(6)、299-301(2017)、 https://www.jstage.jst.go.jp/article/biophys/57/6/57_299/_pdf/-char/ja
*248 Andreas Zamanos et al.、HydraProt: A New Deep Learning Tool for Fast and Accurate Prediction of Water Molecule Positions for Protein Structures、https://pubs.acs.org/doi/pdf/10.1021/acs.jcim.3c01559?ref=article_openPDF
*249 Huynh-Hoa Bui et al.、WATGEN: An algorithm for modeling water networks at protein–protein interfaces、https://onlinelibrary.wiley.com/doi/epdf/10.1002/jcc.20751
*250 Daniel J. Sindhikara et al.、Placevent: An algorithm for prediction of explicit solvent atom distribution—Application to HIV-1 protease and F-ATP synthase、https://onlinelibrary.wiley.com/doi/abs/10.1002/jcc.22984
*251 A. Morozenko et al.、Dowser++, a new method of hydrating protein structures、https://onlinelibrary.wiley.com/doi/epdf/10.1002/prot.25081
*252 https://sbcb.bioch.ox.ac.uk/waterdock/
・・・
(*A-77 Ryan Babbush et al.、The Grand Challenge of Quantum Applications、https://arxiv.org/pdf/2511.09124
)
・・・
(*129 Hakan Doga et al.、A Perspective on Protein Structure Prediction Using Quantum Computers、https://pubs.acs.org/doi/full/10.1021/acs.jctc.4c00067)
(*130 https://newsroom.clevelandclinic.org/2024/06/06/cleveland-clinic-ibm-and-the-hartree-centre-collaborate-to-advance-healthcare-and-life-sciences-through-artificial-intelligence-and-quantum-computing)
・・・
*253 https://meetiqm.com/blog/iqm-resonance-boosts-photodynamic-cancer-therapy/
*254 https://algorithmiq.fi/quantum-ai/
*255 Yanbing Zhou et al.、Quantum Algorithms for Photoreactivity in Cancer-Targeted Photosensitizers、https://arxiv.org/pdf/2512.15889
*256 張元豪、(修士論文)初期状態に依存した量子回路最適化手法と超伝導量子トリットにおける制御ゲート実装手法の開発(改定版)、https://www.icepp.s.u-tokyo.ac.jp/download/master/m2021_jang.pdf
*257 https://qunasys.com/news/posts/quneais_j/
*258 https://www.mitsui.com/jp/ja/topics/2025/1251846_14877.html
*259 https://www.ionq.com/news/ionq-and-ccrm-announce-strategic-quantum-biotech-collaboration-to-accelerate
*260 John Penuel et al.、Detailed assessment of calculating drag force with quantum computers: Explicit time-evolution precludes exponential advantage for nonlinear differential equations、https://arxiv.org/pdf/2406.06323
*261 The Lattice Boltzmann Method、https://link.springer.com/book/10.1007/978-3-319-44649-3
*262 David Jennings et al.、An end-to-end quantum algorithm for nonlinear fluid dynamics with bounded quantum advantage、https://arxiv.org/pdf/2512.03758
*263 Francesco Turro et al.、Practical Application of the Quantum Carleman Lattice Boltzmann Method in Industrial CFD Simulations、https://arxiv.org/pdf/2504.13033
*264 Alexander M. Dalzell、A shortcut to an optimal quantum linear system solver 、https://arxiv.org/pdf/2406.12086
*265 https://www.psiquantum.com/news-import/airbus
*266 David Jennings et al.、Simulating non-trivial incompressible flows with a quantum lattice Boltzmann algorithm、https://arxiv.org/pdf/2512.05781
*267 I.Novikau and I.Joseph、Globalizing the Carleman linear embedding method for nonlinear dynamics、https://arxiv.org/pdf/2510.15715