GitVita
New repositoryImport repositoryNew issue
Your repositories Your dashboard Sign in

microsoft/qdk

Public

mirrored fromhttps://github.com/microsoft/qdkAvailable

Watch0Fork0Star0
CodeCommitsIssuesPull requestsActionsInsightsSecurity
billti/qdk_package

Branches

Tags

  • No tags available.
0Branches0Tags
Go to file
Add file
Code

Clone

HTTPS

Download ZIP

qdk/library/chemistry/src

library/chemistry/src/MixedStatePreparation.qs

308lines · modecode

RawDownload
Latest commit unavailable.
1// Copyright (c) Microsoft Corporation.
2// Licensed under the MIT License.
3
4export MixedStatePreparationRequirements;
5export PurifiedMixedState;
6export PurifiedMixedStateRequirements;
7
8import Std.Arrays.*;
9import Std.Diagnostics.Fact;
10import Std.Math.*;
11import Std.Convert.IntAsDouble;
12import Std.Arithmetic.ApplyIfGreaterLE;
13import Std.StatePreparation.PrepareUniformSuperposition;
14
15import Generators.MultiplexOperationsFromGenerator;
16import Utils.RangeAsIntArray;
17
18/// # Summary
19/// Represents a particular mixed state that can be prepared on an index
20/// and a garbage register.
21///
22/// # Input
23/// ## Requirements
24/// Specifies the size of the qubit registers required to prepare the
25/// mixed state represented by this UDT value.
26/// ## Norm
27/// Specifies the 1-norm of the coefficients used to define this mixed
28/// state.
29/// ## Prepare
30/// An operation that, given an index register, a data register, and a
31/// garbage register initially in the $\ket{0}$, $\let{00\dots 0}$, and
32/// $\ket{00\dots 0}$ states (respectively),
33/// prepares the state represented by this UDT value on those registers.
34struct MixedStatePreparation {
35 Requirements : MixedStatePreparationRequirements,
36 Norm : Double,
37 Prepare : ((Qubit[], Qubit[], Qubit[]) => Unit is Adj + Ctl),
38}
39
40/// # Summary
41/// Represents the number of qubits required in order to prepare a given
42/// mixed state.
43///
44/// # Input
45/// ## NTotalQubits
46/// The total number of qubits required by the represented state preparation
47/// operation.
48/// ## NIndexQubits
49/// The number of qubits required for the index register used by the
50/// represented state preparation operation.
51/// ## NGarbageQubits
52/// The number of qubits required for the garbage register used by the
53/// represented state preparation operation.
54struct MixedStatePreparationRequirements {
55 NumTotalQubits : Int,
56 NumIndexQubits : Int,
57 NumGarbageQubits : Int,
58}
59
60/// # Summary
61/// Returns an operation that prepares a a purification of a given mixed state.
62/// A "purified mixed state" refers to states of the form |ψ⟩ = Σᵢ √𝑝ᵢ |𝑖⟩ |garbageᵢ⟩ specified by a vector of
63/// coefficients {𝑝ᵢ}. States of this form can be reduced to mixed states ρ ≔ 𝑝ᵢ |𝑖⟩⟨𝑖| by tracing over the "garbage"
64/// register (that is, a mixed state that is diagonal in the computational basis).
65///
66/// See https://arxiv.org/pdf/1805.03662.pdf?page=15 for further discussion.
67///
68/// # Description
69/// Uses the Quantum ROM technique to represent a given density matrix,
70/// returning that representation as a state preparation operation.
71///
72/// In particular, given a list of $N$ coefficients $\alpha_j$, this
73/// function returns an operation that uses the Quantum ROM technique to
74/// prepare an approximation
75/// $$
76/// \begin{align}
77/// \tilde\rho = \sum_{j = 0}^{N - 1} p_j \ket{j}\bra{j}
78/// \end{align}
79/// $$
80/// of the mixed state
81/// $$
82/// \begin{align}
83/// \rho = \sum_{j = 0}^{N-1} \frac{|\alpha_j|}{\sum_k |\alpha_k|} \ket{j}\bra{j},
84/// \end{align}
85/// $$
86/// where each $p_j$ is an approximation to the given coefficient $\alpha_j$
87/// such that
88/// $$
89/// \begin{align}
90/// \left| p_j - \frac{ |\alpha_j| }{ \sum_k |\alpha_k| } \right| \le \frac{\epsilon}{N}
91/// \end{align}
92/// $$
93/// for each $j$.
94///
95/// When passed an index register and a register of garbage qubits,
96/// initially in the state $\ket{0} \ket{00\cdots 0}$, the returned operation
97/// prepares both registers into the purification of $\tilde \rho$,
98/// $$
99/// \begin{align}
100/// \sum_{j=0}^{N-1} \sqrt{p_j} \ket{j}\ket{\text{garbage}_j},
101/// \end{align}
102/// $$
103/// such that resetting and deallocating the garbage register enacts the
104/// desired preparation to within the target error $\epsilon$.
105///
106/// # Input
107/// ## targetError
108/// The target error $\epsilon$.
109/// ## coefficients
110/// Array of $N$ coefficients specifying the probability of basis states.
111/// Negative numbers $-\alpha_j$ will be treated as positive $|\alpha_j|$.
112///
113/// # Output
114/// An operation that prepares $\tilde \rho$ as a purification onto a joint
115/// index and garbage register.
116///
117/// # Remarks
118/// The coefficients provided to this operation are normalized following the
119/// 1-norm, such that the coefficients are always considered to describe a
120/// valid categorical probability distribution.
121///
122/// # Example
123/// The following code snippet prepares an purification of the $3$-qubit state
124/// $\rho=\sum_{j=0}^{4}\frac{|\alpha_j|}{\sum_k |\alpha_k|}\ket{j}\bra{j}$, where
125/// $\vec\alpha=(1.0, 2.0, 3.0, 4.0, 5.0)$, and the target error is
126/// $10^{-3}$:
127/// ```qsharp
128/// let coefficients = [1.0, 2.0, 3.0, 4.0, 5.0];
129/// let targetError = 1e-3;
130/// let purifiedState = PurifiedMixedState(targetError, coefficients);
131/// using (indexRegister = Qubit[purifiedState.Requirements.NIndexQubits]) {
132/// using (garbageRegister = Qubit[purifiedState.Requirements.NGarbageQubits]) {
133/// purifiedState.Prepare(LittleEndian(indexRegister), [], garbageRegister);
134/// }
135/// }
136/// ```
137///
138/// # References
139/// - [Encoding Electronic Spectra in Quantum Circuits with Linear T Complexity](https://arxiv.org/abs/1805.03662)
140/// Ryan Babbush, Craig Gidney, Dominic W. Berry, Nathan Wiebe, Jarrod McClean, Alexandru Paler, Austin Fowler, Hartmut Neven
141function PurifiedMixedState(targetError : Double, coefficients : Double[]) : MixedStatePreparation {
142 let nBitsPrecision = -Ceiling(Lg(0.5 * targetError)) + 1;
143 let positiveCoefficients = Mapped(AbsD, coefficients);
144 let (oneNorm, keepCoeff, altIndex) = QuantumROMDiscretization(nBitsPrecision, positiveCoefficients);
145 let nCoeffs = Length(positiveCoefficients);
146 let nBitsIndices = Ceiling(Lg(IntAsDouble(nCoeffs)));
147
148 let op = PrepareQuantumROMState(nBitsPrecision, nCoeffs, nBitsIndices, keepCoeff, altIndex, [], _, _, _);
149 let qubitCounts = PurifiedMixedStateRequirements(targetError, nCoeffs);
150 return new MixedStatePreparation { Requirements = qubitCounts, Norm = oneNorm, Prepare = op };
151}
152
153operation PrepareQuantumROMState(
154 nBitsPrecision : Int,
155 nCoeffs : Int,
156 nBitsIndices : Int,
157 keepCoeff : Int[],
158 altIndex : Int[],
159 data : Bool[][],
160 indexRegister : Qubit[],
161 dataQubits : Qubit[],
162 garbageRegister : Qubit[]
163) : Unit is Adj + Ctl {
164 let garbageIdx0 = nBitsIndices;
165 let garbageIdx1 = garbageIdx0 + nBitsPrecision;
166 let garbageIdx2 = garbageIdx1 + nBitsPrecision;
167 let garbageIdx3 = garbageIdx2 + 1;
168
169 let altIndexRegister = garbageRegister[0..garbageIdx0 - 1];
170 let keepCoeffRegister = garbageRegister[garbageIdx0..garbageIdx1 - 1];
171 let uniformKeepCoeffRegister = garbageRegister[garbageIdx1..garbageIdx2 - 1];
172 let flagQubit = garbageRegister[garbageIdx3 - 1];
173 let dataRegister = dataQubits;
174 let altDataRegister = garbageRegister[garbageIdx3...];
175
176 // Create uniform superposition over index and alt coeff register.
177 PrepareUniformSuperposition(nCoeffs, indexRegister);
178 ApplyToEachCA(H, uniformKeepCoeffRegister);
179
180 // Write bitstrings to altIndex and keepCoeff register.
181 let unitaryGenerator = (nCoeffs, idx -> WriteQuantumROMBitString(idx, keepCoeff, altIndex, data, _, _, _, _));
182 MultiplexOperationsFromGenerator(unitaryGenerator, indexRegister, (keepCoeffRegister, altIndexRegister, dataRegister, altDataRegister));
183
184 // Perform comparison
185 ApplyIfGreaterLE(X, uniformKeepCoeffRegister, keepCoeffRegister, flagQubit);
186
187 let indexRegisterSize = Length(indexRegister);
188
189 // Swap in register based on comparison
190 let register = indexRegister + dataRegister;
191 let altRegister = altIndexRegister + altDataRegister;
192 for i in IndexRange(indexRegister) {
193 Controlled SWAP([flagQubit], (register[i], altRegister[i]));
194 }
195}
196
197// Classical processing
198// This discretizes the coefficients such that
199// |coefficient[i] * oneNorm - discretizedCoefficient[i] * discretizedOneNorm| * nCoeffs <= 2^{1-bitsPrecision}.
200function QuantumROMDiscretization(bitsPrecision : Int, coefficients : Double[]) : (Double, Int[], Int[]) {
201 let oneNorm = PNorm(1.0, coefficients);
202 let nCoefficients = Length(coefficients);
203 Fact(bitsPrecision <= 31, $"Bits of precision {bitsPrecision} unsupported. Max is 31.");
204 Fact(nCoefficients > 1, "Cannot prepare state with less than 2 coefficients.");
205 Fact(oneNorm >= 0.0, "State must have at least one coefficient > 0");
206
207 let barHeight = 2^bitsPrecision - 1;
208
209 mutable altIndex = RangeAsIntArray(0..nCoefficients - 1);
210 mutable keepCoeff = Mapped(
211 coefficient -> Round((AbsD(coefficient) / oneNorm) * IntAsDouble(nCoefficients) * IntAsDouble(barHeight)),
212 coefficients
213 );
214
215 // Calculate difference between number of discretized bars vs. maximum
216 mutable bars = 0;
217 for idxCoeff in IndexRange(keepCoeff) {
218 bars += keepCoeff[idxCoeff] - barHeight;
219 }
220
221 // Uniformly distribute excess bars across coefficients.
222 for idx in 0..AbsI(bars) - 1 {
223 keepCoeff w/= idx <- keepCoeff[idx] + (bars > 0 ? -1 | + 1);
224 }
225
226 mutable barSink = [];
227 mutable barSource = [];
228
229 for idxCoeff in IndexRange(keepCoeff) {
230 if keepCoeff[idxCoeff] > barHeight {
231 barSource += [idxCoeff];
232 } elif keepCoeff[idxCoeff] < barHeight {
233 barSink += [idxCoeff];
234 }
235 }
236
237 for rep in 0..nCoefficients * 10 {
238 if Length(barSink) > 0 and Length(barSource) > 0 {
239 let idxSink = Tail(barSink);
240 let idxSource = Tail(barSource);
241 barSink = Most(barSink);
242 barSource = Most(barSource);
243
244 keepCoeff w/= idxSource <- keepCoeff[idxSource] - barHeight + keepCoeff[idxSink];
245 altIndex w/= idxSink <- idxSource;
246
247 if keepCoeff[idxSource] < barHeight {
248 barSink += [idxSource];
249 } elif keepCoeff[idxSource] > barHeight {
250 barSource += [idxSource];
251 }
252 } elif Length(barSource) > 0 {
253 let idxSource = Tail(barSource);
254 barSource = Most(barSource);
255 keepCoeff w/= idxSource <- barHeight;
256 } else {
257 return (oneNorm, keepCoeff, altIndex);
258 }
259 }
260
261 return (oneNorm, keepCoeff, altIndex);
262}
263
264/// # Summary
265/// Returns the total number of qubits that must be allocated
266/// in order to apply the operation returned by
267/// `PurifiedMixedState`.
268///
269/// # Input
270/// ## targetError
271/// The target error $\epsilon$.
272/// ## nCoefficients
273/// The number of coefficients to be specified in preparing a mixed state.
274///
275/// # Output
276/// A description of how many qubits are required in total, and for each of
277/// the index and garbage registers used by the
278/// `PurifiedMixedState` function.
279function PurifiedMixedStateRequirements(targetError : Double, nCoefficients : Int) : MixedStatePreparationRequirements {
280 Fact(targetError > 0.0, "targetError must be positive");
281 Fact(nCoefficients > 0, "nCoefficients must be positive");
282
283 let nBitsPrecision = -Ceiling(Lg(0.5 * targetError)) + 1;
284 let nIndexQubits = Ceiling(Lg(IntAsDouble(nCoefficients)));
285 let nGarbageQubits = nIndexQubits + 2 * nBitsPrecision + 1;
286 let nTotal = nGarbageQubits + nIndexQubits;
287 return new MixedStatePreparationRequirements { NumTotalQubits = nTotal, NumIndexQubits = nIndexQubits, NumGarbageQubits = nGarbageQubits };
288}
289
290operation WriteQuantumROMBitString(idx : Int, keepCoeff : Int[], altIndex : Int[], data : Bool[][], keepCoeffRegister : Qubit[], altIndexRegister : Qubit[], dataRegister : Qubit[], altDataRegister : Qubit[]) : Unit is Adj + Ctl {
291 if keepCoeff[idx] >= 0 {
292 ApplyXorInPlace(keepCoeff[idx], keepCoeffRegister);
293 }
294 ApplyXorInPlace(altIndex[idx], altIndexRegister);
295 if Length(dataRegister) > 0 {
296 for i in IndexRange(data[idx]) {
297 if data[idx][i] {
298 X(dataRegister[i]);
299 }
300 }
301 for i in IndexRange(data[altIndex[idx]]) {
302 if data[altIndex[idx]][i] {
303 X(altDataRegister[i]);
304 }
305 }
306 }
307}
308
309