量子コンピュータで自由に初期状態を作る方法

量子コンピュータの基礎:振幅の初期化 - yukobaのブログの続き。

量子コンピュータでは、各状態の振幅の絶対値が確率な訳であり、観測するとどれかが決まるので、つまり、確率ベクトルです。どうやって確率ベクトルの初期状態を作るかは大事です。上の日記は、(1/3, 1/3, 0, 1/3) みたいに、全て同じ確率の物しか作れなかったので、(4/7, 2/7, 0, 1/7)のように、ばらばらの確率で作れるように自分で改良しました。つまり、\sqrt{\frac{4}{7}}|0> + \sqrt{\frac{2}{7}}|1> + \sqrt{\frac{1}{7}}|3>の作り方です。当然のことながら、確率の合計は1にする必要があります。これに加えて、組み込み関数のCPhase()を使うと、自由に偏角も設定できます。

qcl -i initAmpGeneral.qcl で、testInitAmp() で実行できます。実行には、QCL - A Programming Language for Quantum Computersが必要です。

: SPECTRUM x: <0,1>
0.57143 |0>, 0.28571 |1>, 0.14286 |3>

このような結果が返ってきます。4/7 = 0.57143 です。

ソースコード

変わったのは、initAmpStateGen() とそれを呼び出している部分です。

// Initializing amplitude distribution (振幅分布の初期化)
// 元の論文は確率 1 / #patterns にしかできなかったのを、任意の確率で初期化できるようにしたもの
// 論文の上位ビットと下位ビットの順番が一般の順序と逆であることに注意!特にcのビット順が問題。

// 論文の3〜6行目
qufunct initAmpFlip(qureg x, qureg c, int pattern, int prevPattern) {
  int j;
  // 論文の5,6行目はF^0であり、CNotはF^1なので、ビット反転。
  Not(c[0]);
  // 論文の3〜5行目
  for j = 0 to #x - 1 {
    if bit(pattern, j) xor bit(prevPattern, j) { CNot(x[j], c[0]); }
  }
  // 論文の6行目
  CNot(c[1], c[0]);
  // 反転したc[0]を元に戻す
  Not(c[0]);
}

// 論文の7行目
operator initAmpStateGen(qureg c, real p) {
  real p1; real p2;
  // p1^2 + p2^2 = 1 より行列はユニタリ
  p1 = sqrt(1 - p);
  p2 = sqrt(p);

  Matrix4x4(
    1, 0, 0, 0,
    0, 1, 0, 0,
    0, 0, p1, -p2,
    0, 0, p2, p1,
    c);
}

// c1 == i1 かつ c2 == i2 の時に、a を反転する。i1, i2 = {0, 1}。#a, #c1, #c2 = 1。
qufunct fredkin(qureg a, qureg c1, boolean i1, qureg c2, boolean i2) {
  if not i1 { Not(c1); }
  if not i2 { Not(c2); }
  CNot(a, c1 & c2);
  if not i2 { Not(c2); }
  if not i1 { Not(c1); }
}

// 論文の8〜10行目
// x が z と等しければ、g[#x - 1] のフラグを立てる
qufunct initAmpAnd(qureg x, quvoid g, int z) {
  int k;
  fredkin(g[0], x[0], bit(z, 0), x[1], bit(z, 1));
  for k = 2 to #x - 1 {
    fredkin(g[k - 1], x[k], bit(z, k), g[k - 2], true);
  }
}

// 振幅分布の初期化
operator initAmp(quvoid x, int vector patterns, real vector values) {
  int i; int prevPattern;
  real leftProb;
  qureg g[#x - 1];
  qureg c[2];
  
  prevPattern = 0;
  leftProb = 1;
  for i = 0 to #patterns - 1 {
    // 論文の3〜6行目
    initAmpFlip(x, c, patterns[i], prevPattern);
    prevPattern = patterns[i];
    
    // S^p。論文の7行目。
    initAmpStateGen(c, values[i] / leftProb);
    leftProb = leftProb * (1 - values[i] / leftProb);
    
    // save。パターンに等しいやつだけの c[1]をフリップ
    // 論文の8〜10行目(AND)
    initAmpAnd(x, g, patterns[i]);
    // 論文の11行目
    CNot(c[1], g[#x - 2]);
    // 論文の12〜14行目(AND逆戻し)
    !initAmpAnd(x, g, patterns[i]);
  }
  // 論文の15行目
  Not(c[0]);
}

procedure testInitAmp() {
  const patternBitSize = 2;
  int vector patterns = vector(0, 1, 3);
  real vector values = vector(4.0/7.0, 2.0/7.0, 1.0/7.0);
  qureg x[patternBitSize];
  initAmp(x, patterns, values);
  dump x;
}