Биномиальное распределение с параметрами n и p в теории вероятностей — распределение количества «успехов» в последовательности из n независимых случайных экспериментов, таких, что вероятность «успеха» в каждом из них постоянна и равна p.

Рассмотрим случай, когда p=0.5 - это делается только для упрощения примера.

В этом случае, согласно теории, вероятность выпадения исхода k на вероятностном пространстве из целых чисел равно P(k)=C(k,n)/2**n, где C(k,n) = n!/(k!(n-k)!) - коэффициент бинома Ньютона.

Поставим перед собой цель - сформировать в массиве кубитов, который мы будем рассматривать как регистр из нескольких кубитов, состояние SUM SQRT(C(k,n))|k>

Из теории вероятности имеем, что если x1,...,xn - конечная последовательность независимых случайных величин, и имеющих одинаковое распределение Бернулли с параметром p, то есть при каждом i=1..n величина 1 («успех») и 0 («неудача») с вероятностями p и q=1-p соответственно, тогда случайная величина y=x1+...+xn имеет биноминальное распределение с параметрами n и p.

В квантовых вычислениях, мы имеем факт, что применение трансформации Адамара H к кубиту в состоянии |0> даёт нам его в равновероятном состоянии для исходов |0> и |1>, то есть в состоянии |0>+|1>

Воспользуемся этими фактами, и сформируем значение нужного нам регистра, как "сумму значений" n независимых кубитов, находящихся в состоянии |0>+|1>

Таким образом алгоритм прост:

1. Берём несколько кубитов в состоянии |0> (регистр), для хранения значения в диапазоне [0..n] в двоичном представлении,

2. Берём n кубитов xi в состоянии |0>

3. Применяем преобразование Адамара к xi

4. Суммируем xi с накоплением суммы в регистре (да, нам потребуется реализовать аналог двоичной арифметики, но на кубитах)

5. Очищаем и освобождаем кубиты xi

В результате данных преобразований, получаем регистр в состоянии SUM SQRT(C(k,n))|k>

Для получения конкретного значения с биноминальным распределением, нам надо коллапсировать регистр (то есть провести измерение значений его кубитов)
и преобразовать полученный вектор из |0> и |1> в целое число (двоичное представление числа)

Перейдём к кодированию с помощью QDK от Microsoft

1. Воспользуемся инструкциями с сайта microsoft.com по установки QDK под VS Code (https://learn.microsoft.com/ru-ru/training/modules/qsharp-create-first-quantum-development-kit/2-install-quantum-development-kit-code).

2. Создадим новый проект (или клонируйте мой https://github.com/dprotopopov/qbinom)

3. Добавим

   3.1. Трансформацию "инкремента регистра"

   /// Реализация операции инкремента, то есть трансформации |k> -> |k+1>
   
       operation Inc(register: Qubit[]) : Unit is Ctl {
           let n = Length(register);
           for idx in 1..n {
               Controlled X(register[0..n-idx-1], register[n-idx]);
           }
       }

   3.2. Трансформацию "сумма битов"

   /// Реализация операции суммирования, то есть трансформации |abcde...> -> |a+b+c+d+e+...>
   
       /// Это не самая эффективная реализация, но самая простая для кодирования
   
       operation Sum(qubits: Qubit[], register: Qubit[]) : Unit {
           for q in qubits {
               Controlled Inc([q], register);
           }
       } 

   3.3. И, соответственно, операцию заполнения регистра суммой значений случайных битов

   /// Реализация операции биноминального распределения, то есть n -> SUM SQRT(C(k,n))|k>
   
       operation Binom(n: Int, register: Qubit[]) : Unit {
           use qubits = Qubit[n] {
               ApplyToEach(H, qubits);
               Sum(qubits, register);
               ResetAll(qubits);
           }
       }

4. Создадим пример для проверки

   4.1. то есть несколько раз повторим эксперимент

   4.1.1. Заполним регистр из кубитов состоянием SUM SQRT(C(k,n))|k>

   4.1.2. Коллапсируем регистр, и получим конкретное значение k

   4.1.3. Запишем, что нам выпал исход k

   4.2. Выведем статистику выпадения исходов

        // Параметры нашего примера

        let n = 10;
        let tests = 1000;
        let k = BitSizeI(n);

        use register = Qubit[k] {

            // Аллокируем массив для подсчёта количества исходов экспериментов

            mutable arr = [0, size = n + 1];

            // Повторям эксперимент несколько раз, с подсчётом колличества полученных исходов

            for _ in 1..tests {

                // Установим в register состояние SUM|k>

                Binom(n, register); 

                // Для получения конкретного значения необходимо измерить значения кубиртов в регистре

                // и преобразовать полученный результат (вектор из |0> или |1>) в целое  число

                let results = ForEach(M, register);
                let i = ResultArrayAsInt(results);

                // Добавим полученное значение в счётчик

                set arr w/= i <- arr[i] + 1;

                // Очищаем кубиты после использования

                ResetAll(register);
            }

            // Выводим полученную таблицу частот 

            for s in arr {
                let p = 100 * s / tests;
                Message($"{p}%");
            }
        }
}

namespace qbinom {

    open Microsoft.Quantum.Canon;
    open Microsoft.Quantum.Intrinsic;
    open Microsoft.Quantum.Arrays;
    open Microsoft.Quantum.Math;
    open Microsoft.Quantum.Measurement;
    open Microsoft.Quantum.Convert;
    
    /// Реализация операции инкремента, то есть трансформации |k> -> |k+1>
    operation Inc(register: Qubit[]) : Unit is Ctl {
        let n = Length(register);
        for idx in 1..n {
            Controlled X(register[0..n-idx-1], register[n-idx]);
        }
    }

    /// Реализация операции суммирования, то есть трансформации |abcde...> -> |a+b+c+d+e+...>
    /// Это не самая эффективная реализация, но самая простая для кодирования
    operation Sum(qubits: Qubit[], register: Qubit[]) : Unit {
        for q in qubits {
            Controlled Inc([q], register);
        }
    }

    /// Реализация операции биноминального распределения, то есть n -> SUM SQRT(C(k,n))|k>
    operation Binom(n: Int, register: Qubit[]) : Unit {
        use qubits = Qubit[n] {
            ApplyToEach(H, qubits);
            Sum(qubits, register);
            ResetAll(qubits);
        }
    }

    @EntryPoint()
    operation Main() : Unit {
        Message("Hello quantum world!");

        // Параметры нашего примера
        let n = 10;
        let tests = 1000;


        let k = BitSizeI(n);
        use register = Qubit[k] {

            // Аллокируем массив для подсчёта количества исходов экспериментов
            mutable arr = [0, size = n + 1];

            // Повторям эксперимент несколько раз, с подсчётом колличества полученных исходов
            for _ in 1..tests {
                // Установим в register состояние SUM SQRT(C(k,n))|k>
                Binom(n, register); 

                // Для получения конкретного значения необходимо измерить значения кубиртов в регистре
                // и преобразовать полученный результат (вектор из |0> или |1>) в целое  число
                let results = ForEach(M, register);
                let i = ResultArrayAsInt(results);

                // Добавим полученное значение в счётчик
                set arr w/= i <- arr[i] + 1;

                // Очищаем кубиты после использования
                ResetAll(register);
            }

            // Выводим полученную таблицу частот 
            for s in arr {
                let p = 100 * s / tests;
                Message($"{p}%");
            }
        }
    }
}

Итог

PS C:\Projects\qbinom> dotnet run
Hello quantum world!
0%
1%
4%
11%
20%
25%
20%
11%
3%
1%
0%

Похоже? Мне кажется, что - да.

Ссылки

https://github.com/dprotopopov/qbinom