microsoft/qdk

Public

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

CodeCommitsIssuesPull requestsActionsInsightsSecurity
7c4c1b4c424591b52bf0c1009540c321c7756e69

Branches

Tags

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

Clone

HTTPS

Download ZIP

samples/estimation/EkeraHastadFactoring.qs

195lines · modecode

1/// # Sample
2/// Resource Estimation for Integer Factoring
3///
4/// # Description
5/// In this sample we concentrate on costing quantum part in the algorithm for
6/// factoring RSA integers based on Ekerå and Håstad
7/// [ia.cr/2017/077](https://eprint.iacr.org/2017/077) based on the
8/// implementation described in
9/// [arXiv:1905.09749](https://arxiv.org/abs/1905.09749). This makes it ideal
10/// for use with the Azure Quantum Resource Estimator.
11import Std.Convert.*;
12import Std.Math.*;
13import Std.ResourceEstimation.*;
14import Std.Arrays.*;
15import Std.Arithmetic.*;
16import Std.TableLookup.*;
17
18// !!! IMPORTANT !!!
19// When computing resource estimtes from the VS Code plugin directly on this
20// file, make sure that you set the error budget to 0.333.
21
22operation Main() : Unit {
23 // Try different instances of the algorithm by commenting in and out
24 // the following lines. You can find more RSA numbers at
25 // https://en.wikipedia.org/wiki/RSA_numbers
26
27 // RSA-100 (330 bits)
28 EkeraHastad(330, 1522605027922533360535618378132637429718068114961380688657908494580122963258952897654000350692006139L, 7L);
29
30 // RSA-1024 (1024 bits)
31 // EkeraHastad(1024, 135066410865995223349603216278805969938881475605667027524485143851526510604859533833940287150571909441798207282164471551373680419703964191743046496589274256239341020864383202110372958725762358509643110564073501508187510676594629205563685529475213500852879416377328533906109750544334999811150056977236890927563L, 7L);
32
33 // RSA-2048 (2048 bits)
34 // EkeraHastad(2048, 25195908475657893494027183240048398571429282126204032027777137836043662020707595556264018525880784406918290641249515082189298559149176184502808489120072844992687392807287776735971418347270261896375014971824691165077613379859095700097330459748808428401797429100642458691817195118746121515172654632282216869987549182422433637259085141865462043576798423387184774447920739934236584823824281198163815010674810451660377306056201619676256133844143603833904414952634432190114657544454178424020924616515723350778707749817125772467962926386356373289912154831438167899885040445364023527381951378636564391212010397122822120720357L, 7L);
35}
36
37/// # Summary
38/// Main algorithm based on quantum phase estimation
39///
40/// # Reference
41/// [ia.cr/2017/077, Section 4.3](https://eprint.iacr.org/2017/077)
42operation EkeraHastad(numBits : Int, N : BigInt, g : BigInt) : Unit {
43 let x = ExpModL(g, ((N - 1L) / 2L), N);
44 let xinv = InverseModL(x, N);
45
46 let m = numBits / 2;
47 use c1 = Qubit[2 * m];
48 use c2 = Qubit[m];
49 use target = Qubit[numBits];
50
51 let ne = 3 * m;
52 let cpad = Ceiling(2.0 * Lg(IntAsDouble(numBits)) + Lg(IntAsDouble(ne)) + 10.0);
53
54 // The algorithm uses the coset representation to replace modular
55 // addition inside MultiplyExpMod with regular addition.
56 use padding = Qubit[cpad];
57
58 ApplyToEach(H, c1 + c2);
59 InitCoset(N, target, padding);
60 MultiplyExpMod(g, N, c1, target + padding);
61 MultiplyExpMod(xinv, N, c2, target + padding);
62
63 Adjoint ApplyQFT(c1);
64 Adjoint ApplyQFT(c2);
65}
66
67// ------------------------------ //
68// Modular arithmetic (constants) //
69// ------------------------------ //
70
71/// Window size for exponentiation (c_exp)
72internal function ExponentWindowLength_() : Int { 5 }
73
74/// Window size for multiplication (c_mul)
75internal function MultiplicationWindowLength_() : Int { 5 }
76
77// ------------------------------- //
78// Modular arithmetic (operations) //
79// ------------------------------- //
80
81/// # Summary
82/// Encodes register in coset representation
83///
84/// # Reference
85/// [arXiv:quant-ph/0601097, Section 4.1](https://arxiv.org/abs/quant-ph/0601097)
86internal operation InitCoset(mod : BigInt, xs : Qubit[], padding : Qubit[]) : Unit {
87 use helper = Qubit();
88 let cpad = Length(padding);
89 let n = Length(xs);
90
91 let combined = xs + padding;
92
93 for j in 0..cpad - 1 {
94 Controlled IncByLUsingIncByLE([helper], (RippleCarryCGIncByLE, mod, combined[j..j + n - 1]));
95
96 ApplyIfLessOrEqualL(X, mod, combined[j..j + n - 1], helper);
97 }
98}
99
100/// # Summary
101/// Computes zs *= (base ^ xs) % mod (for a large register xs)
102///
103/// # Reference
104/// [arXiv:1905.07682, Fig. 7](https://arxiv.org/abs/1905.07682)
105internal operation MultiplyExpMod(
106 base : BigInt,
107 mod : BigInt,
108 xs : Qubit[],
109 zs : Qubit[]
110) : Unit {
111 let expWindows = Chunks(ExponentWindowLength_(), xs);
112
113 within {
114 RepeatEstimates(Length(expWindows));
115 } apply {
116 let i = 0; // in simulation this i must be iterated over IndexRange(expWindows)
117
118 let adjustedBase = ExpModL(base, 1L <<< (i * ExponentWindowLength_()), mod);
119 MultiplyExpModWindowed(adjustedBase, mod, expWindows[i], zs);
120 }
121}
122
123/// # Summary
124/// Computes zs *= (base ^ xs) % mod (for a small register xs)
125///
126/// # Reference
127/// [arXiv:1905.07682, Fig. 6](https://arxiv.org/abs/1905.07682)
128internal operation MultiplyExpModWindowed(
129 base : BigInt,
130 mod : BigInt,
131 xs : Qubit[],
132 zs : Qubit[]
133) : Unit {
134 let n = Length(zs);
135
136 use qs = Qubit[n];
137 AddExpModWindowed(base, mod, 1, xs, zs, qs);
138 AddExpModWindowed(InverseModL(base, mod), mod, -1, xs, qs, zs);
139 for i in IndexRange(zs) {
140 SWAP(zs[i], qs[i]);
141 }
142}
143
144/// # Summary
145/// Computes zs += ys * (base ^ xs) % mod (for small registers xs and ys)
146///
147/// # Reference
148/// [arXiv:1905.07682, Fig. 5](https://arxiv.org/abs/1905.07682)
149///
150/// # Remark
151/// Unlike in the reference, this implementation uses regular addition
152/// instead of modular addition because the target register is encoded
153/// using the coset representation.
154internal operation AddExpModWindowed(
155 base : BigInt,
156 mod : BigInt,
157 sign : Int,
158 xs : Qubit[],
159 ys : Qubit[],
160 zs : Qubit[]
161) : Unit {
162 // split factor into parts
163 let factorWindows = Chunks(MultiplicationWindowLength_(), ys);
164
165 within {
166 RepeatEstimates(Length(factorWindows));
167 } apply {
168 let i = 0; // in simulation this i must be iterated over IndexRange(factorWindows)
169
170 // compute data for table lookup
171 let factorValue = ExpModL(2L, IntAsBigInt(i * MultiplicationWindowLength_()), mod);
172 let data = LookupData(factorValue, Length(xs), Length(factorWindows[i]), base, mod, sign, Length(zs));
173
174 use output = Qubit[Length(data[0])];
175
176 within {
177 Select(data, xs + factorWindows[i], output);
178 } apply {
179 RippleCarryCGIncByLE(output, zs);
180 }
181 }
182}
183
184internal function LookupData(factor : BigInt, expLength : Int, mulLength : Int, base : BigInt, mod : BigInt, sign : Int, numBits : Int) : Bool[][] {
185 mutable data = [[false, size = numBits], size = 2^(expLength + mulLength)];
186 for b in 0..2^mulLength - 1 {
187 for a in 0..2^expLength - 1 {
188 let idx = b * 2^expLength + a;
189 let value = ModulusL(factor * IntAsBigInt(b) * IntAsBigInt(sign) * (base^a), mod);
190 set data w/= idx <- BigIntAsBoolArray(value, numBits);
191 }
192 }
193
194 data
195}
196