Abstract
The security of RSA is generally understood to be based on the complexity of factoring, and key size parameters are extrapolated from the general number field sieve (GNFS). However, this may not accurately represent RSA security in practical scenarios.
An under-appreciated 2007 algorithm of Joux, Naccache, and Thomé allows an attacker to forge RSA signatures after temporary access to a raw RSA signing/decryption oracle in time close to the special number field sieve (SNFS) without factoring the key. We implement and run this algorithm for 1024-bit RSA. In total, the attack took 1380 CPU core-years over five calendar months, and made 2³² oracle queries. Most of this time is precomputation; after the precomputation the attacker can forge any signature of choice, offline, in 180 core-years.
We carried out our attack using a hardware security module (HSM) as the signing oracle, thus demonstrating the ability to impersonate the HSM through black-box API interactions, without exfiltrating the key. Blind RSA schemes also provide such a signing oracle.
Extrapolating our empirical running times to larger key sizes, we conclude that the concrete security of RSA with a signing oracle should be 15 to 30 bits lower than the factoring-based security estimates for the 1024-bit to 4096-bit RSA parameters that are common in practice. Even 4096-bit RSA does not appear to meet a 128-bit security level in this attack model. This highlights a gap in current RSA-type security assumptions, and gives classical cryptanalytic evidence in favor of moving away from RSA entirely during the current post-quantum transition.
1 Introduction
In this paper, we carry out an attack that allows us to forge arbitrary 1024-bit RSA digital signatures (or decrypt arbitrary ciphertexts) using academic computational resources, without factoring the modulus. This attack reduces RSA security levels below acceptable thresholds for modern cryptographic deployments, and demonstrates a theoretical gap in RSA-based security assumptions and practical parameter choices.
Alternative title: Nearly Special-Number-Field-Sieve-Speed Signature Forgery Sans Factoring N (NSNFSSSFSFN)
The threat model we consider is a “lunchtime attack,” or a non-adaptive chosen ciphertext attack (IND-CCA1)—an attacker who obtains temporary access to an RSA signature/decryption oracle for a given key, and then is challenged to generate a signature forgery or ciphertext decryption after losing access. We demonstrate that temporary oracle access, together with much more modest computational resources than factoring the modulus would require, can give the attacker the permanent ability to forge arbitrary signatures for that key pair offline. That is, the attacker can effectively use oracle queries to steal abilities equivalent to possessing the private key, without actually computing it.
The attack runs in two separate stages: an initial precomputation and oracle query phase, followed by a fully offline signature forgery/decryption phase that can be repeated at will for arbitrary targets and requires no further oracle queries. The precomputation stage depends only on the modulus N and took us about 1,200 core-years to complete. The oracle query phase required us to query around 2³² raw, unpadded RSA signatures/decryptions. We carried out our oracle queries using an HSM (Hardware Security Module). The final offline phase takes around 180 core-years to forge any signature (or target decryption). In comparison, current estimates for factoring a 1024-bit RSA modulus would take 500,000–1,000,000 CPU core-years with current implementations [15]. The computation time for our attack can almost certainly be improved, and should be interpreted as an upper bound.
The attack we implement is from a 2007 paper of Joux, Naccache, and Thomé [42]. Their attack runs in asymptotic time (L_N(1/3, 1.577 + o(1))), close to that of the special number field sieve. This is significantly less than the (L_N(1/3, 1.923 + o(1))) running time of the general number field sieve for factoring. We implemented several algorithmic improvements and had to overcome a number of computational challenges to carry out this attack for 1024-bit RSA. We detail the implications for both practical and theoretical settings below.
1.1 Practical impacts
RSA digital signatures are still in wide use for important protocols like TLS and OAuth [13], with 2048-bit keys the most common choice. RSA’s days are numbered, however: NIST plans to deprecate it by 2030 and disallow it by 2035 in favor of post-quantum algorithms [54], and the transition to hybrid and post-quantum algorithms is underway. Although NIST disallowed 1024-bit RSA in 2013 for signature generation [6] and modern implementations have largely migrated to 2048 bits, 1024-bit keys are still allowed for signature validation under legacy use [7], and remain in use today. (See Section 7.4.)
While the hardness of factoring is generally used to determine RSA key strengths (e.g. in NIST SP 800-56B Rev. 2 [53]), there are several practical usage scenarios that match the signing oracle access attack model we consider in this paper. Such an attack likely does not pose an immediate operational threat to most modern RSA deployments, but it demonstrates that RSA provides significantly less security than the factoring-based security estimates would suggest in these settings.
Security estimates. 1024-bit RSA is generally considered to provide 80 bits of security; we estimate that the algorithm we implement has 2⁶⁵ running time for 1024-bit RSA. Scaling the asymptotic running time to 2048-bit RSA gives a 2⁹⁰ running time for the attack, with 2⁴³ queries, a significant reduction from the 112-bit security that 2048-bit RSA is generally understood to provide. Even 4096-bit RSA fails to achieve 128 bits of security in this attack model—we extrapolate that such an attack would take 2¹¹⁹ computation time and require 2⁵⁷ queries. See Section 6 for more details.
Hardware security modules. Hardware security modules (HSMs) are often used to protect secret keys and accelerate cryptographic operations in enterprise environments. Current security auditing best practices, as well as standards such as PCI DSS [57], encourage or require the use of HSMs to protect key material. Standards such as FIPS 140-3, ISO/IEC 19790, ISO/IEC 24759, and ISO/IEC 15408 (Common Criteria) define requirements for these types of cryptographic modules to be tamper-proof and include side-channel countermeasures to protect keys even against an attacker with physical access.
The PKCS 11 standard [55] implemented by nearly all HSMs provides an interface for raw RSA operations that exactly matches the oracle query functionality required by our attack. Thus, in a threat model requiring an HSM to protect private keys in such a scenario, our attack can be viewed as allowing an attacker to effectively extract capabilities equivalent to the private key in a feasible amount of time with only black-box access to the HSM itself.
Blind signature schemes. RSA-based blind signature schemes also provide the type of signature query oracle required by our attack, and thus this attack should be interpreted as lowering the bit security of RSA in protocols using blind signatures. We analyze some deployed applications in Section 7.2 and observe that the 2048-bit key sizes common in practice provide too low of a security threshold.
1.2 Theoretical impacts
The problem of breaking RSA (that is, computing eth roots modulo an integer N of unknown factorization) has not been shown to be equivalent to factoring N in the nearly 50 years since the invention of the RSA cryptosystem. Results in various models [1] have not been shown to apply to (Z^*_N).
There are a number of computational assumptions associated with RSA security. A prominent one for signatures is the one-more RSA assumption [8].
However, we observe that the attack we implement, which we name (\tilde{\sqrt{\mathrm{NFS}}}), does not align with one-more RSA or variants despite being a clear practical threat.
We propose that there is a gap both in the practice of using the running time of GNFS to set parameter sizes for these RSA-type assumptions, and in the assumptions themselves. We do not resolve this gap. We suggest that precomputation attacks should be captured in RSA-type assumptions, and that subexponential-time attacks (including their lower-order terms) should be reasonably analyzable in those frameworks. We have named our attack setting “delayed-target RSA” in line with the naming scheme of Koblitz and Menezes [45]. We discuss this in detail in Section 5.
1.3 RSA progress in general
The algorithmic result we use has been known for two decades but the impact on RSA parameter choices does not appear to have been understood in the community. We hope that our concrete implementation for a landmark key size makes the difference much clearer. 1024-bit RSA has long been considered a minimal security level, and remains in use today in some specialized cases.
Coincidentally, the factorization of RSA-260, along with that of a 1277-bit special modulus, was publicly posted [48] a few days after our own computation finished. At least two individuals had been independently working on RSA-260, using AI to translate CADO-NFS to run on GPUs, and then using GPU compute from their AI company employers to run it [70]. Both our findings and the pre-2026 factoring records were completed using academic computational resources and engineering skills, a far cry from industry- or state-level adversaries.
That leaves substantial leeway before one can get a sense of the feasibility limit that would apply to adversaries with larger computing power. Extrapolations that RSA-1024 could be factored in weeks using industry-grade GPU facilities are no fantasy. The important takeaway is that improvements to sieving for the GNFS, including to ECM or sieving on GPUs, can definitely be obtained. Any speed improvements to CADO-NFS would apply in a straightforward way to the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm as well. Weis and Lu report that factoring RSA-260 cost 2⁵ to 13.5 GPU-years (respectively) using their implementations. A back-of-the-envelope scaling suggests that our 1024-bit RSA signature forgery computation would take between 2 and 9 GPU-years if we had used their implementations; for 2048 bits a very rough estimate would be 2²⁷ GPU-years.
NFS development. CADO-NFS [64] is an open-source implementation of the number field sieve algorithms for factoring and discrete logarithm that has been used for several recent computational records. The algorithm that we implement did not have a public implementation prior to our work. Our implementation makes use of the highly optimized CADO-NFS code for sieving, linear algebra, and descent, among other algorithmic procedures, but uses them in a different setting, in a different order and with different parameters.
Our engineering work does not hugely transform the CADO-NFS codebase we started with—at least for the usual number field sieve algorithms. However it took significant technical effort to port CADO-NFS to our setting, and our work could be useful for other such extensions to the NFS. The progression of our 1024-bit computation—from parameter selection to CADO-NFS roadblocks to algorithmic changes—will be an informative data point for future large computations. In particular, we implemented interesting variations of algorithms for descent, and for eth roots in a number field.
The future of RSA. The community has been focusing on quantum computing attacks as a reason to retire RSA, and the classical cryptanalysis of RSA has been largely regarded as settled for decades, with occasional factoring records set by teams that track with general-purpose computing improvements over the years. We hope that our computation here illustrates that there are still surprises to be found in the classical cryptanalysis of RSA, and if more attention had been paid there would have been even stronger reasons to retire RSA earlier.
The security of RSA and the family of NFS-based algorithms used to calculate RSA key sizes should not be viewed as a static monolith. Even the classical security of RSA is more complex than is generally appreciated in the community.
While we do not believe that the attack we implement should be cause for immediate alarm among practitioners managing cryptographic inventory in industry—it is practical in an academic sense rather than in the script-kiddie sense, and the attack model of a raw, unpadded RSA signing oracle is very strong in practice even if academically relevant—we believe that this attack demonstrates that 2048 to 4096-bit RSA falls short of modern cryptographic security levels requiring 128 bits of security or larger. Given this, the current transition to post-quantum algorithms presents an opportunity to transition to modern algorithms and security levels in the next few years.
1.4 Summary of computational results
We summarize the results of our large computations in Figure 1. Our code is available at https://github.com/ucsd-hacc/NSNFSSSFSFN.
| Key length |
Precomputation (core-years) |
Queries |
Individual forgery (core-years) |
| 512 |
0.01 |
2²³ |
0.001 |
| 666 |
3.5 |
2²⁷ |
0.07 |
| 768 |
56 |
2²⁷ |
42 |
| 1024 |
1200 |
2³² |
180 |
Figure 1. Empirical running times for (\tilde{\sqrt{\mathrm{NFS}}}). We experimented with different parameter settings giving different tradeoffs for different phases of the algorithm.
| Algorithm |
Complexity |
Modulus size |
Year, citation |
CPU time |
Wall time |
| GNFS-DL |
G |
795 |
2019 [14] |
3200 core-years |
|
| GNFS |
G |
795 |
2019 [14] |
1000 core-years |
|
| GNFS |
G |
829 |
2020 [14] |
2700 core-years |
|
| SNFS-DL |
S |
1024 |
2016 [28] |
385 core-years |
2 months |
| SNFS |
S |
1039 |
2007 [3] |
225 core-years |
> 6 months |
| GNFS (GPU) |
G |
862 |
2026 [48] |
13.5 GPU-years |
|
| GNFS (GPU) |
G |
896 |
2026 [71] |
TBA |
|
| SNFS (GPU) |
S |
1277 |
2026 [48] |
6.6 GPU-years |
|
| (\tilde{\sqrt{\mathrm{NFS}}}) |
S+ |
512 |
2007 [42] |
> 2 core-days |
|
| (\tilde{\sqrt{\mathrm{NFS}}}) |
S+ |
512 |
this work |
4.8 core-days |
70 minutes |
| (\tilde{\sqrt{\mathrm{NFS}}}) |
S+ |
666 |
this work |
3.6 core-years |
33 hours |
| (\tilde{\sqrt{\mathrm{NFS}}}) |
S+ |
768 |
this work |
98 core-years |
3 weeks |
| (\tilde{\sqrt{\mathrm{NFS}}}) |
S+ |
1024 |
this work |
1380 core-years |
5 months |
Figure 2. Related record computations for factoring, discrete logarithm, and (\tilde{\sqrt{\mathrm{NFS}}}). We categorize complexity as either General (G), Special (S), or slightly above special (S+).
2 Preliminaries
Number field sieve. The fastest known classical general-purpose algorithms for factoring and discrete log are a family of Number Field Sieve (NFS) algorithms. To disambiguate, we will refer specifically to the GNFS and SNFS variants, and their discrete log variants GNFS-DL and SNFS-DL. The “G” variants work for a general modulus, meaning any generic integer modulus. The “S” variants only work for a special modulus—in brief, a modulus that admits an efficiently discoverable number field where sieving is faster.
The topic of this paper is what we call the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm (“eth root NFS”). It was introduced in 2007 by Joux, Naccache, and Thomé (JNT) [42, §4]. (\tilde{\sqrt{\mathrm{NFS}}}) is not a straightforward extension of the usual GNFS or GNFS-DL algorithms. While it reuses large parts of them, the steps are slightly different, and the parameter tradeoffs differ even further. We will elaborate on this later; for now we stress that (\tilde{\sqrt{\mathrm{NFS}}}) should be viewed as a different algorithm in the NFS family. For a readable but exhaustive resource on the NFS, we refer the reader to the lecture notes and references in [69].
Figure 2 lists some record computations along with our own results. An up-to-date list is often kept on Wikipedia [72]. We have omitted the breakdown of precomputation/individual computation in this table; there is a similar division for the discrete log variants but not for factoring.
Subexponential complexity. The running times of all of the algorithms we discuss in this paper are subexponential, with running time expressed asymptotically as (L_N(\alpha, c)^{1+o(1)}) with
[
L_N(\alpha, c) = \exp(c \cdot (\ln N)^\alpha(\ln\ln N)^{1-\alpha}).
]
In this paper, we will always encounter (\alpha = 1/3) with a varying (c).
The running time of the GNFS for factoring is (L_N(1/3, (64/9)^{1/3} + o(1))) with N the integer to be factored, and (c \approx 1.923). This is the fastest known asymptotic running time to factor generic integers, such as well-generated RSA moduli. The running time of the GNFS-DL for prime field discrete logarithm is the same, with N the prime modulus.
The running time of the SNFS for factoring is (L_N(1/3, (32/9)^{1/3} + o(1))), which improves (c \approx 1.526). The SNFS can be used to factor integers that are efficiently expressable in a special form; for example it applies to integers close to a power of 2. To give an intuitive sense of this difference, a 1024-bit RSA modulus has yet to be factored publicly, but a 1039-bit Mersenne number was publicly factored using the SNFS algorithm in 2007 with 225 core-years of compute effort [3]. A similar gap exists in the discrete logarithm case [28].
For our purposes, “nearly special” complexity is (L_N(1/3, 1.577 + o(1))). This is the overall complexity of the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm when a practical variant called extension sieving is used. (Which we did.) The full statement of the running time, along with the proof, were given in [42].
Our problem setting. The most apt name we have seen for our RSA problem of interest is “delayed-target RSA” following the terminology of Koblitz and Menezes [45] when discussing JNT’s work. We give some more detail and context for this in Section 5. For now, we describe the problem as follows. An attacker is assumed to operate in two phases, sharing state between them. In the first phase, they are given an RSA public key (N, e) and may query an oracle which, on input w, returns (w^{1/e} \bmod N). (In usual notation, (1/e = d \bmod \phi(N)).) In the second phase, the attacker receives a target (t \leftarrow$ Z^*_N) and loses access to the oracle. They must produce (t^{1/e} \bmod N) to win the game.
It turns out that, as currently known, an attack on delayed-target RSA does not yield an attack for usual RSA (without the oracle), nor for factoring. There is a more subtle relationship to one-more RSA [8], a more popular assumption used to prove security of blind signatures. We discuss that relationship in Section 5.
Other similar attack settings. Problem definitions like delayed-target RSA have made some appearances in the literature. On the attack side, JNT considered our delayed-target RSA as well as variants for other padding schemes [42]. Joux, Lercier, Naccache and Thomé [41] considered delayed-target CDH in finite-field groups, where a similar attack applies. Koblitz and Menezes [45] gave attacks for delayed-target DL on hyperelliptic curves. On the constructive side, delayed-target CDH was used to prove security of Freeman’s identification scheme [27]. Delayed-target DDH has appeared in proofs for the CCA1 security of ElGamal [47] and Cramer-Shoup [9]. We are not aware of any proofs relying on delayed-target RSA; if any exist they are likely under a different name.
3 The (\tilde{\sqrt{\mathrm{NFS}}}) algorithm, with our additions
The (\tilde{\sqrt{\mathrm{NFS}}}) algorithm is from Section 4 of JNT [42]. That paper refers to the “one-more RSA inversion” problem but to be technically precise, it is attacking delayed-target RSA as we defined in the Preliminaries.
Our implementation of (\tilde{\sqrt{\mathrm{NFS}}}) uses CADO-NFS [64], SageMath [67], PARI [66] and additional code written in Python. It makes heavy use of existing NFS ideas and folklore. The major modifications and new techniques are discussed in this section, along with enough background to be understandable. The specific parameters and progression of the 1024-bit computation are discussed in Section 4.
3.1 Bird’s-eye view
The (\tilde{\sqrt{\mathrm{NFS}}}) algorithm can be viewed as a malleability attack, combining ideas from the GNFS and GNFS-DL algorithms.
Broadly, the goal is to be able to express any target t as a product of elements whose eth roots we know from precomputation and oracle queries. Of course, if these “elements” are only rational primes, the running time and success probability are wholly impractical. Rather we develop an equivalence, following the same notations as in [42, §4.1]:
[
t^d =
\left(\frac{u}{v}S(m)\right)^d R(m)\cdot T(m)^d\cdot U(m)^d
]
where m is from our setup; u, v, and S are from descent; T is from precomputation; U is from linear algebra; R is from taking a root in the number field; and all the dth powers are known from precomputed oracle queries. Roughly, lowercase letters are integers and uppercase letters are rational functions.
In a nutshell, the algorithm approaches SNFS complexity because it only needs to sieve for smooth relations on one side, compared to two sides for GNFS.
The “other side” is traded for oracle queries. However, the descent phase does not become easier, which necessitates rebalancing the steps of the computation.
3.2 Detailed view: Precomputation and oracle queries
Polynomial selection. Given a modulus N, we select two polynomials f, g such that (f(m) = g(m) = 0 \bmod N) for a shared root m. For all relevant cases we ran, an optimal choice of degrees is (\deg(g) = 1) and (\deg(f) \in {5, 6}). Our number field is then (K = \mathbb{Q}(\alpha)) where (f(\alpha) = 0). Its ring of integers is (O_K). For later steps it is useful to select f such that the public RSA exponent e is unramified (or even inert) in K, and there are several other small inert primes. Other than these extra checks, the usual NFS polynomial selection can be used.
There are assumed to be three smoothness parameters input to the algorithm: (B_{\mathrm{alg}}, B_{\mathrm{ext}},) and (B_{\mathrm{rat}}). Let:
[
F_{\mathrm{alg}} = {p \mid p \subset O_K \text{ is prime and } |p| < B_{\mathrm{alg}}}
]
[
F_{\mathrm{ext}} = {p \mid p \subset O_K \text{ is prime and } B_{\mathrm{alg}} < |p| < B_{\mathrm{ext}}}
]
[
F_{\mathrm{rat}} = {p \mid p \in \mathbb{Z} \text{ is prime and } p < B_{\mathrm{rat}}}.
]
We call these sets “factor bases.”
One-sided precomputation sieving. We made slight fixes to the CADO-NFS las program to allow for one-sided sieving. For (\tilde{\sqrt{\mathrm{NFS}}}), in contrast to the GNFS and the SNFS, precomputation sieving is done only on the algebraic side. It is broken down into two steps. In the first, algebraic sieving, sieving is performed to find relations that are (B_{\mathrm{alg}})-smooth. By “relation” we mean a tuple (a, b) and a factorization such that
[
(a-b\alpha) = \prod p_i^{\mu_i}
]
where ((a-b\alpha)) is the ideal generated by (a-b\alpha \in K), and the (p_i) are all in (F_{\mathrm{alg}}). We find enough of these algebraic relations such that every (p_i) in (F_{\mathrm{alg}}) is included in at least one relation.
The second sieving step is extension sieving. For every (q \in F_{\mathrm{ext}}), special-q sieving is performed to find one relation where the factorization includes exclusively q and primes (p_i) that are in the smaller (F_{\mathrm{alg}}) (or at least of norm less than the norm of q). In practice, this does not need to occur for every (q \in F_{\mathrm{ext}}), but it should occur for a good portion of them.
Oracle queries. From the above precomputation, we make algebraic queries for all of the algebraic relations, and extension queries for all of the extension relations. That is, for a relation (a, b) we make a query to learn ((a-bm)^d \in Z^*_N).
We also make rational queries to learn (p_i^d \in Z^*N) for all rational primes (p_i) from 2 to (B{\mathrm{rat}}). If our polynomial g is not monic (which is nearly always the case), we make queries for all the prime factors of its leading coefficient. Often, these are already covered by the rational queries.
3.3 Detailed view: Descent
Three-stage descent. Upon receiving a forgery target t, the goal of descent is to produce an initial split (u/v = t \bmod N), and a rational function S, such that:
[
(S(\alpha)) = \prod_{p_i \in F_{\mathrm{alg}}\cup F_{\mathrm{ext}}} p_i^{\mu_i}
]
and
[
u/v \cdot S(m) = \prod_{p_j \in F_{\mathrm{rat}}} p_j^{\nu_j}.
]
In other words, the ideal generated by (S(\alpha)) is smooth over the entire algebraic side, and (u/v \cdot S(m)) is smooth over the rational side. The smoothness bounds in these products vary as the descent tree grows, but they will reach (B_{\mathrm{ext}}) and (B_{\mathrm{rat}}) by the last step of descent.
On the surface, we reused sieving code for descent that the GNFS-DL would use in CADO-NFS. However a difference arises in parameter tradeoffs. In the GNFS-DL and SNFS-DL, the same smoothness bounds are used in descent and in precomputation. Given the massive amount of precomputation, the expected time to find a two-sided relation must remain quite low. The (\tilde{\sqrt{\mathrm{NFS}}}) setting is different. Descent must find a few hundred two-sided relations, but one-sided sieving over its millions of relations can get away with a much more aggressive smoothness bound. This tradeoff was not discussed by JNT.
For a concrete number, in our 1024-bit computation, we selected (B_{\mathrm{ext}} = 2^{35}) and (B_{\mathrm{rat}} = 2^{36}). Each two-sided descent relation took several hours to find, while each one-sided precomputation relation took less than a second. We could have sped up descent by increasing (B_{\mathrm{ext}}) to (2^{36}), but this would have doubled the amount of extension sieving. Descent took less than one month, while extension sieving took 2.5 months. The cost of doubling extension sieving was frightening.
(i) Descent initialization. Our descent algorithm takes place in three stages. The initialization stage borrows ideas from the 1024-bit SNFS-DL computation in 2016 [28]. It includes a las-plus-ECM step and a “large q” step. In the former step, candidate initial splits (u/v = t \bmod N) are found by sieving for candidate factorizations of u and v, while allowing for large (we chose 365-bit) cofactors. Those cofactors are handled by a more efficient ECM setup than is built into las. We aimed for 135-bit smoothness here. The “large q” step handles the rational primes of size about 95 to 135 bits. They are too large to input to las as a special-q, so we used a technique that we call “embedded-polynomial sieving.” We describe it below in step (iii), where we used the same technique.
(ii) Descent middle. The middle of descent is its most tractable stage. We were able to use the existing las descent code in CADO-NFS. For all of our computations, we required that this stage begin with special-q values of at most 95 bits. This is because the sieving lattice in las uses 64-bit integers, and has some dependence on the skewness of the polynomial.
On its own, this simple implementation of the descent middle was enough to handle moduli up to a 666-bit computation. However, we needed the additional initialization and bottom stages for the 768 and 1024-bit computations.
(iii) Descent bottom. To implement the bottom of descent, we used two strategies. The first is composite sieving and the second is embedded-polynomial sieving. During our large computations, we found that while both eventually find relations, composite sieving was usually faster. Composite sieving is simple to describe: rather than perform a single special-q (or q) job, it is possible to perform in parallel many special-((qq')) jobs for small primes (or composites) (q'). The benefit is parallelism. Concretely, we would run one las command, but with a long list of special-((qq')) inputs on which to perform multi-threaded sieving. Composite sieving (albeit in a different form) was used in the computation from [14], which used CADO-NFS. It took very little work to adapt it to our setting.
Embedded-polynomial sieving is an idea from the “large q” step and that we repurposed from [28]. Rather than find a special-q (or q) relation (a, b) under polynomials f, g it is possible find any relation ((a', b')) under transformed polynomials (f_q, g_q) such that ((a', b')) can be turned into a relation (a, b) under f, g which includes the prime q as a factor. The process is a homographic transformation. We briefly detail how that works, choosing as an example a special-q on the algebraic side: choose ((a_0, b_0)) and ((a_1, b_1)) two basis vectors of the lattice of (a, b)’s such that q divides ((a-b\alpha)). Here, the lattice is to be reduced according to the skewness parameter of the polynomial f. Let now
[
f_q = f\left(\frac{a_0x+a_1}{b_0x+b_1}\right)(b_0x+b_1)^\delta
]
(likewise for (g_q), with (\delta) the relevant degree). Notice that by construction, the skewness of (f_q) becomes approximately 1 after this transformation. The benefits here are parallelism—since many small primes can be used to sieve under (f_q, g_q)—as well as injecting variety into the sieving lattices, in case a particular q has odd behavior. Although these ideas are not built into the las siever, these bootstrapped (but folklore) uses of las were critical to making descent faster.
3.4 Detailed view: Linear algebra
The linear algebra step takes the rational function S output by descent and the information from precomputation, and it builds a matrix-vector equation where the solution yields rational functions U, T such that (S(\alpha)T(\alpha)U(\alpha) \in K) is an eth power. That eth power is the input to the final root stage. While most of this step was handled by the CADO-NFS program bwc and extensive parameter optimization, we made two other nontrivial changes.
Linear algebra optimizations. In contrast to factoring, where the matrix is mod 2, and to discrete log, where the matrix is mod (\ell) for a large factor (\ell) of the group order, our implementation constructs a matrix that is mod the RSA public exponent e. We chose (e = 65537) to match common use. Our matrix equation is inhomogeneous, like discrete log, and unlike factoring. These differences led to extensive debugging of existing CADO-NFS linear algebra code. One interesting optimization we were able to use involves large matrix multiplications. The modulus (e = 65537) is small enough that we could represent entries as float64 values, and use optimized sparse matrix multiplication libraries [36, 65, 61] while keeping enough precision to reduce mod e at the end. This worked well up to our 768-bit computation, but the 1024-bit computation required more precision.
Filtering additions. Filtering is technically part of precomputation, but it involves the matrix we just discussed above. In this matrix, columns are prime ideals in (F_{\mathrm{alg}}), rows are algebraic relations, and entries are valuations mod 65537. CADO-NFS includes a step called filtering which performs elementary operations on the matrix to make a smaller linear system. In this filtered matrix, the rows are actually linear combinations of the original relations, and the columns are a small subset of (F_{\mathrm{alg}}). We used “killer relations” to map between filtered row or column vectors, to vectors over the full pre-filtered space. This is similar to a step called “reconstruction” in GNFS-DL but we needed our own implementation for this use case. The basic idea is to iteratively cancel ideals in a prescribed order, at each step using a so-called killer relation that includes only the current ideal and previously cancelled ideals. During our precomputation filtering phase, we thus built up this ordering of killer relations, and used it to appropriately apply or remove filtering on the linear algebra vectors.
3.5 Detailed view: eth root in a number field
Instead of using a single algorithm to find an eth root in a number field, we took a three-stage approach. This is the final step of (\tilde{\sqrt{\mathrm{NFS}}}) before assembling the final forgery using Eq. (1).
(i) Parallelized Montgomery root. On the 768 and 1024-bit computations, we found that a coarse but more parallel version of Montgomery’s algorithm worked well for the first 90% of norm reduction. We followed the approach described in [68], but rather than selecting generators by analyzing the log-embeddings, we simply chose the first generator of the desired norm. We could choose up to a few dozen generators in parallel, and run only one accumulation step on the product of all the generators. Experimentally, this worked well in our computations. Moreover, it did not work well to try to pick generators intelligently (using the approach described below) in independent parallel instances; the log-embedding improvements tended to be cancelled out or even reversed. We believe there could be a more interesting way to exploit parallelism in Montgomery’s algorithm. This root algorithm can also be used (mod 2) in GNFS for factoring.
(ii) Sequential Montgomery root. This matches the approach described in [68]. Here, we ran a single thread to perform one reduction step at a time, selecting generators to balance the log-embeddings as much as possible. A purely sequential Montgomery root was fast for everything up to our 666-bit computation, but it was outperformed by the three-stage approach for 768 and 1024 bits.
(iii) Hybrid CRT-p-adic root. The Montgomery root, given a field element (R(\alpha)^e), produces a list of (\gamma_i) such that
[
R(\alpha)^e = \delta^e \prod_i \gamma_i^e
]
for some unknown (\delta) with bounded coefficients; with each Montgomery iteration, the bound on the coefficients of (\delta) ideally becomes smaller. One hopes that eventually (\delta) will become 1, and (R(\alpha) = \prod_i \gamma_i), but in practice this is not always the case. We then have to deal with the leftover (\delta^e), which we get as a product expression
[
\prod_i(a_i-b_i\alpha)^{m_i}\prod_j\gamma_j^{-e}.
]
The product is much too large to fully evaluate, but it can be evaluated mod some modulus p, as this prevents the coefficients from blowing up.
For the 768-bit computation and below, we handled the leftover case using the p-adic lifting approach described in [68]. In short, we choose an inert prime p (i.e. f is irreducible mod p) so that (F_p[x]/(f(x)) \simeq F_{p^{\deg f}}), a finite field where eth roots can be taken easily. We evaluate the product expression for (\delta^e) mod (p^k) (where (p^k) is above the coefficient bound on (\delta), but far below that of (\delta^e)), take its eth root mod p (since eth roots in (F_{p^{\deg f}}) are easy), and then Hensel lift to recover the coefficients of (\delta \bmod p^k). Since (p^k) is above the coefficient bound, this recovers (\delta).
For our 1024-bit computation, we parallelized this step by using many inert primes (p_i) and solving for the coefficients of (\delta \bmod p_i^{k_i}) in parallel. We then recover (\delta) via CRT reconstruction. (Here (k_i) can be smaller because we only need (\prod_i p_i^{k_i}) to be above the bound, not each (p_i^{k_i}).) This is thus a hybrid of the p-adic lifting and the CRT approaches of [68]. Compared to a purely-CRT approach, the hybrid approach let us use fewer primes (since we only require (\prod_i p_i^{k_i}) be above the bound, not (\prod_i p_i)). As pointed out by Bernard, Fouque, and Lesavourey [12], unlike in the square root case, for odd e there is no need to guess the correct signs of the eth roots mod each (p_i) before CRT reconstruction: eth roots are unique in (F_{p_i^{\deg f_i}}) as long as e is relatively prime to (p_i^{\deg f_i}-1).
The Montgomery root algorithm is asymptotically faster than CRT and/or p-adic, but the latter is much more parallelizable. It seems likely that the root step could be made faster by carefully balancing the two.
3.6 Detailed view: Parameter selection
To the frustration of the reader and any interested replicators, our parameters were chosen through a combination of intuition and small-scale experiments. For many parameters, it was infeasible to exhaustively run experiments on every combination of allowed parameters on academic resources. Our process was to first decide on high-level choices ((B_{\mathrm{rat}}, B_{\mathrm{ext}}, B_{\mathrm{alg}})) then run experiments with mid-level choices (A, mfb, etc.). For sieving, our experiments tended to optimize the number of unique relations found per second. For descent sieving, we also accounted for RAM usage; we increased A, mfb until parallelized sieving ran out of memory on our 512GB machines. For the Montgomery root, our experiments increased the per-iteration norm bound until we either ran out of precision in our called SageMath and PARI code, or ran out of memory.
Regarding the smoothness bounds, (B_{\mathrm{alg}}) was always selected first, and not changed afterwards. For some of the computations, we ran descent and extension sieving on different subsets of machines, so that we could compare their expected finishing times and adjust as needed. In particular, descent sieves from top (large primes) to bottom (smaller primes). We ran extension sieving starting from (B_{\mathrm{alg}}) (small primes) and going up to (B_{\mathrm{ext}}) (larger primes). Thus both bounds (B_{\mathrm{ext}}, B_{\mathrm{rat}}) could be increased without losing past work.
We ran fewer experiments for linear algebra but were informed by past NFS computations. We chose to follow [28] which selected a matrix of primes bounded by (2^{31}) for a 1024-bit SNFS-DL computation. This was successfully run for a much larger matrix modulus, so we expected to slightly outperform their timing. We ended up taking 19 days; theirs took one month.
4 1024-bit computation details
The public RSA key on our HSM consisted of:
e = 65537
N = aa8bb728b43ba939dc282ae5eace9e939ae5a985398db5816d4e266a1a64986d
21a19152fb3cb738aa7ecd3b868f5ed42d6afd2ce33d450ab68ff0d8705604e4
72d8a235a06b743bb5ea2b731585cafe38d9b5ec6f0f1aefd86dc26426891144
67bb7ede0a3314daed0f3d184a7f605105ebd49dee447bbc23e331dfa8e1c113
with (\log_2(N) = 1024). The random target that was selected for forgery, and the eth root that we ultimately computed, were:
t = 804fbc7eb701fdc8fa1147900eb44395ace53ac96a9c676726131ef4fb4158fa
33c50e3bf4bda303d8a6a8c51c0a3031b6a0c12378d7eb3b760ad760c30790af
d64517920348c12e26ed75f233a4e731f772309c5df0873d8bcab2c116e3ff04
3ea1ee3379eac327c58686b40f87edae312b0534a50df24206caa1e220368eb8
t1/e = 944678c4efa56b54dcaff74fac3cca0acc5dd571d18af653171693d8560fc2cf
752e62edc5ba2a6f378e7afdebf3d3ffb3af6e1c139764d3753e9976da2dad12
aa94e8365ae7cb599fdc528bc0618fc789fc0a27aac8403e101d4611c5faff62
5832ca3c3509a7f021553d49e434dade3a097b0adfed94468a1d5e1d9ee68f88
The target was selected at random from (Z^*_N). It can be easily verified that ((t^{1/e})^e \bmod N = t).
Our computation ran from March 13 to August 31, 2026. Due to breaks and (un)expected cluster downtime, the actual running time was about 5 months. We chose to run certain steps in parallel to more quickly identify problems that arise only for very large moduli; there were such problems at nearly every step. If running only one step at a time on 100% of the cluster, the breakdown of the 5 months would be about: 2.5 months of extension sieving, 0.5 months of algebraic sieving, 1 month of descent, and 1 month of linear algebra and the final steps.
| Step |
CPU time (h) |
Wall time (h) |
Notes |
| Precomputation |
1200 core-years |
|
|
| Polynomial selection |
39,992.01 |
20.17 |
|
| Algebraic sieving |
1,311,975.08 |
2,760.00 |
Small portion of cluster |
| Filtering* |
24.08 |
7.85 |
|
| Extension sieving |
9,240,757.82 |
2,640.00 |
Ran on > 1/2 of cluster |
| HSM queries |
4.1 billion queries |
|
|
| Rational queries* |
N/A |
123.12 |
2,874,398,515 queries |
| Algebraic queries* |
N/A |
6.59 |
144,973,282 queries |
| Extension queries* |
N/A |
45.15 |
1,048,047,534 queries |
| Individual forgery |
180 core-years |
|
|
| Descent (i) |
17,141.54 |
19.82 |
Ran on < 1/2 of cluster; each core used significant RAM |
| Descent (ii) |
18,811.00 |
268.36 |
|
| Descent (iii) |
1,385,242.92 |
739.21 |
|
| Linear algebra* |
146,263.96 |
447.04 |
|
| Root* |
2,210.93 |
147.16 |
|
Figure 4. Timing of our 1024-bit computation. All steps other than those marked with a () can be trivially parallelized over more machines, meaning the wall-clock time can be easily decreased with more resources.*
4.1 Timing, memory usage, and queries
The timing of each step is summarized in Figure 4. The CPU time should be viewed as additive. The wall-clock time is not, due to how we scheduled steps to be run on different subsets of the cluster. Some steps of the algorithm are trivially parallelizable, meaning that if twice as many machines were used, the wall-clock time would be cut in half. This applies to all steps except those marked with a (*) in Figure 4; notably, linear algebra and the Montgomery root.
Our “cluster” included a few different types of machines. We used:
- Cluster A: 45 machines with 48–88 cores and 256GB–1TB of RAM
- Cluster B: 16 machines with 256 cores and 1TB of RAM
- Cluster C: (Varying) machines with 32+ cores and 192GB+ of RAM
- Cluster D: 2 machines with 512 cores and 1TB+ of RAM
The above machines have hyperthreading enabled and we report virtual cores; there are half as many physical cores. The bulk of the work was done on Cluster A which has relatively dated processors from 2016 but we had uninterrupted access to all machines. We had priority access to Cluster B and some access to Cluster C, where we mostly performed extension sieving. We had short-term access to Cluster D where we performed some extension sieving, as well as high-RAM steps including the Montgomery root.
4.2 Parameters selected
The global smoothness bounds we selected were:
[
B_{\mathrm{rat}} = 2^{36}, \qquad B_{\mathrm{ext}} = 2^{35}, \qquad B_{\mathrm{alg}} = 2^{31}.
]
These are the three most consequential to the running time. For our polynomials, we selected (\deg(g) = 1) as always for a rational side, and (\deg(f) = 6).
Sieving occurs in three different steps in the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm, and they all used different sieving parameters. Ordered from easiest to hardest to manage, they are: algebraic sieving, extension sieving, then descent sieving.
The linear algebra and root parameters described below were determined largely from the above, as well as by the exponent being (e = 65537).
By this point, the reader may notice that we counterintuitively set (B_{\mathrm{rat}} > B_{\mathrm{ext}}). In usual NFS computations with a rational side, one would prefer to set (B_{\mathrm{rat}} < B_{\mathrm{ext}}). Here, recall that (B_{\mathrm{ext}}) is the smoothness bound on the algebraic side, when two-sided sieving is performed. For example, in the 2019 RSA-240 factorization via GNFS [14], (B_{\mathrm{rat}} = 2^{36}) was used, and the equivalent of (B_{\mathrm{ext}}) was set to (B_{\mathrm{ext}} = 2^{37}). Roughly, this is because the algebraic side is “heavier” meaning that for a given relation, the algebraic side tends to have a larger norm, which is less likely to factor into small primes. We did observe the same behavior in two-sided sieving. However, the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm has a different global tradeoff between (B_{\mathrm{rat}}) and (B_{\mathrm{ext}}). To make two-sided sieving in descent tractable, one can either increase (B_{\mathrm{ext}}), which roughly doubles the cost of extension sieving in precomputation and doubles the number of extension queries; or one can increase (B_{\mathrm{rat}}), which incurs no extra sieving cost but still doubles the number of rational queries. We made the second choice more often than the first, because we were more concerned with sieving time than with query time.
4.3 Possible tradeoffs and alternative regimes
We chose the above parameters with the goal of making the end-to-end computation work in a handful of months. Different high-level parameters could be selected to optimize other goals. The point of this discussion is not to provide definitive answers about alternative parameter regimes (which would need to actually be run) but to suggest that ours is not the only one possible.
For example, another realistic goal is to minimize the number of HSM queries, especially if the attacker has limited online access to the HSM. In this paradigm, one could feasibly select (B_{\mathrm{rat}} = 2^{34}) and (B_{\mathrm{ext}} = 2^{34}). Our experiments suggest that this would cut the number of queries to 1.5 billion while significantly increasing the time spent in descent (perhaps by a month).
On the other extreme, a goal could be to minimize the time spent in individual forgery. We believe it would be hard to decrease the time spent in linear algebra, but descent could be sped up by increasing (B_{\mathrm{rat}}) and (B_{\mathrm{ext}}). Our experiments suggest that the individual forgery could be made to take 5–6 weeks, but this starts to require tens of billions of queries.
4.4 Trip report: Precomputation
Polynomial selection. We used a similar polynomial selection approach to prior factorization computations. After using CADO-NFS scripts to generate candidates and perform root optimization, we narrowed our search to the 100 best polynomials based on the “Murphy E” value.
Given these 100 candidates, we performed about 8 hours of sieving experiments on our cluster, testing the polynomials with different sieving parameters. Ranking by the average time to find a unique relation, there was a clear winner and runner-up. We also verified that these polynomials satisfied necessary properties for our root computation to work: first that (e = 65537) was unramified in the corresponding number field; and second that there were several small inert primes with which we could perform p-adic lifting.
A last consideration was skewness. Because descent involves sieving where the special-q or special-q can be up to 95 bits, we needed to be careful with the size of entries in the sieving lattice, which in CADO-NFS are represented as 64-bit unsigned integers. Since the larger entries are of size about (\sqrt q \cdot) skewness, we chose to aim for a skewness of about (2^{16}). This was reflected in our input parameters to the CADO-NFS candidate-generation program. The polynomials we selected are:
g(x) = + 583567398878296587761x
− 583508354875246158635574786673190632633220109293
f(x) = − 3823002681315339817451520x6
− 7350594038360643424557871006920x5
− 216020850394775040895128014268735128x4
+ 1766821832922784779107212285802811346027985x3
− 10216745863874085463847140035876784402919224735x2
− 29161682848128620601357534837213030324970598132169770x
− 1140093999763033664163418590780932917435250519974808516832
Here, g is the rational-side polynomial, f is the algebraic-side polynomial, and (g(m) = f(m) = 0 \bmod N) for a shared root m. The relevant number field is (K = \mathbb{Q}(\alpha)) where (f(\alpha) = 0). The skewness is (2^{67.925}).
Algebraic sieving. Algebraic sieving was the most straightforward step of our computation. Given our polynomial f, we used CADO-NFS to create a factor base of all prime ideals of (O_K) with norm at most (B_{\mathrm{alg}} = 2^{31}). Then we used las to sieve over large ranges of special-q until we had at least one relation covering every q in the factor base.
We used a sieving area of (A = 2^9) and increased to (A = 2^{10}) as we needed more relations. We used a cofactor of mfb = 100 and increased to mfb = 120; in other words, we allowed for 3–4 large primes. We actually used a smaller factor base limit of (\lim = 2^{27}) since we found that this smaller limit (with a larger mfb) was faster on our machines. CADO-NFS allows for multi-layer bucket sieving, which we did use (see [14, §4.2]). We set bkthresh1 = (2^{26}). Recall that this is all one-sided sieving, which explains the perhaps smaller-than-expected parameters.
Filtering (which we detail below) is a step that processes the found relations into a matrix. After filtering, it turned out that some of our factor-base ideals became unlinked, due to singleton elimination or other simplifications. For these ideals, we performed special-q sieving to find a few more relations for each q. Rather than increasing A and mfb, we found that it was more helpful to perform composite sieving over (q \cdot q') for small prime (q'). While some unlinked ideals may not be a problem for the GNFS, for our setting, we wanted to cover 100% of ideals so that any descent output could give us a solvable linear algebra equation.
Extension sieving. Although extension sieving was by far the most computationally expensive step in our computation, it was luckily predictable and consistent to run. Recall that here, we are sieving on one side with special-q that are between (B_{\mathrm{alg}} = 2^{31}) and (B_{\mathrm{ext}} = 2^{35}), with smoothness bound (B_{\mathrm{alg}} = 2^{31}). We used the same parameters as above for algebraic sieving, but started out at the higher end of (A = 2^{10}) and mfb = 120.
Other than the scale and sizes, there were two key differences in extension sieving compared to algebraic sieving. First, the goal is to find only one relation per target ideal. Although we set the CADO-NFS parameter to “exit early” we also needed to de-duplicate afterwards, since this flag does not fully stop all threads. It is a different de-duplication step than the one used for algebraic sieving.
Second, we did not need nor want 100% coverage of all the ideals between (B_{\mathrm{alg}} = 2^{31}) and (B_{\mathrm{ext}} = 2^{35}). On average, relations are straightforward to find, but the tail end of hard-to-include ideals quickly becomes an expensive task. We aimed for 90% coverage. Then, descent jobs would be likely to include our linked ideals, and for descent relations involving the 10% of unlinked ideals, it was only a little more work to find a second relation avoiding these. The coverage of ideals that we ultimately found was all prime ideals whose norm was below (2^{31}), 90% below (2^{33}), and down to 60% below (2^{35}). Practically, this struck a reasonable balance with descent. Overall, however, extension sieving is a unique trick to the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm, and we believe other tradeoffs could be justified.
Filtering and matrix preparation. The final step of precomputation is to prepare the matrix for linear algebra. We used standard CADO-NFS steps, plus an additional step to account for what we called “killer relations.” We used all of the CADO-NFS filtering programs (dup1, dup2, purge, merge-dl, replay-dl) with discrete-log flags, since our linear algebra setup closely matches that of GNFS-DL rather than GNFS. The only nontrivial parameters we set were to aim for a target density of (2^{10}), and to keep at least 100 extra rows after purge. The sparse part of our filtered matrix had 10,557,419 rows and 10,557,319 columns, among which there were 2,217,211,872 nonzero entries. Each row had an additional 6 dense righthand columns for the Schirokauer maps.
4.5 Trip report: Oracle queries
To instantiate the RSA signature oracle, we first had to obtain access to an HSM. Enterprise-grade HSMs are expensive proprietary cryptographic accelerators that are opaque by design. They are usually sold exclusively on the business-to-business (B2B) market and not to end consumers. The smaller, consumer-grade HSMs that we considered typically disallowed 1024-bit RSA. We gained access to the following two models, through different means:
- Thales Luna K6 (PCIe 7000): This HSM is an expansion card that can be plugged into a server’s PCIe slot. It was donated to us and had an estimated market price in the high four digits or low five digits before being deprecated in June 2021 (for companies) and December 2023 (for governments) [50].
- Thales Luna S750: A third party allowed us to use their HSM remotely. It is still sold for $22,000–32,000 US dollars [35, 40] and does not appear to have an end-of-life date yet.
In addition to the above HSMs, we unsuccessfully attempted to set up several other HSMs locally but were unable to find a copy of the proprietary software that is required to set them up. In the following, we first describe how we configured our local Luna K6 HSM, followed by our system architecture for querying both the local and remote HSM. To conclude, we report the wall-clock time for querying, as well as the signature speed of the two HSMs.
Setting up the Luna K6 HSM. After installing the Luna K6 HSM in a slot that is (backwards) compatible with PCIe x4, we had to install the necessary software to interact with the HSM. Together with the Luna K6 HSM, we received a CD-ROM with proprietary software for this HSM model. However, due to changes to the Linux kernel headers and other dependencies, we discovered that Ubuntu 10.04 is the latest OS for which this software can be compiled. We could not emulate this old Ubuntu version since the HSM registers as a PCIe bridge, which cannot be passed through to a VM. Hence, we did a bare metal install, which meant we had to manually backport several packages, for example to support the network card of our server which was produced after Ubuntu 10.04 was declared end-of-life. To interface with the HSM we installed the pycryptoki [62] PKCS#11 wrapper for Python developed by Thales.
For our local HSM, we had to obtain a PIN Entry Device (PED) [60] and PED keys [59] from the second-hand market. The former is an external peripheral device with a small display, numeric keypad, and a few action keys. It serves as a trusted input/output device which bypasses the untrusted server where the HSM is plugged into a PCIe slot. The PED is used during setup time to program PED keys—which are programmable devices of a similar form factor as USB keys—with key material for various privileged roles. We further enabled “user activation” and “auto-activation”, which generate a secret value displayed on the PED that can later be used to start user sessions on the HSM from remote.
In the global HSM settings, we disabled FIPS 140-2 approved operation mode to enable raw RSA operations (referred to as CKM RSA X 509 in the Cryptographic Token Interface Standard [20]). Querying unpadded messages is necessary for our attack. Afterwards, we create a partition which inherits the global HSM settings. We chose to do key generation outside of the HSM and then import it because HSMs wipe their state on a power outage (of which we had multiple during this project). For this, we generate a 1024-bit RSA key in PKCS#1 format and store it in DER encoding and imported it with the Certificate Management Utility (CMU) from the Luna client software. HSM deployments would instead likely use the cloning or wrapping features, through which the keys can be transferred directly to another HSM of the same model without the keys ever leaving the HSM unencrypted.
Query architecture. Our goal was to have a joint system that runs queries both on our Luna K6 HSM and the remote Luna S750. For the latter, we were not given access to the third party’s network but were only allowed to pass software to them. Hence, we designed the following three-piece system architecture.
First, a local query handler groups queries into batches and uploads them to the web server. It also regularly attempts to download new signature batches and aggregates them. Second, the web server offers an authenticated REST API to upload and download batches, metadata, performance metrics, and debug information. Third, the remote query handler fetches query batches from the web server, runs them on the HSM, writes the resulting signatures into batches, and uploads them.
Various code handles network and hardware failures, performs rate-limiting, and re-synchronizes the handlers. Both query handlers send measurements of their upload and download bandwidth to the web server. Additionally, the remote query handler uploads the measured signature signing speed, as well as CPU and memory usage to the web server. All components are highly parallelized. We carefully chose the batch sizes, number of batches in transit, number of threads to write, download, pre or post-process batches, and webserver parallelism such that the bottleneck was at the HSM, not the network bandwidth or batch pipeline. The remote query handler catches most exceptions, reports them to the web server to enable remote debugging, and then restarts gracefully.
Query timing. The HSM does not use parallelism within sessions, even under a high query load. Instead of multiple OS threads using the same session (which throws exceptions on the HSM), a thread per session should be used to benefit from parallelism and maximize the signature speed. The official administration guide recommends using 20–40 threads [30]; we achieved an average of 6500 signatures per second with 40 threads. Computing the first batch of rational queries on our Luna K6 HSM took over a month due to many pauses for debugging corner cases that only occurred at scale (e.g. unlikely concurrency bugs).
Based on the average signature speed, running the rational queries would have taken 123.12 hours without interruptions. The algebraic and extension queries took 6.59 and 45.15 hours respectively, showing that the HSM ran at impressively consistent signature speed. The overhead from our query architecture was overall minor and stemmed from ramping up the pipeline utilization, error handling and re-synchronization in case a garbled batch arrived, and post-processing for the very last batches.
Additionally, we ran the same queries on the remote Luna S750 HSM, which was located on a different continent. We initially expected the Luna S750 HSM to be faster than ours, since it is a newer generation, but that was not the case. The rational, algebraic, and extension queries took 317.25, 7.08, and 118.93 hours, respectively. This corresponds to a signature speed of approximately 2500 sig/s for the rational and extension queries, while the algebraic queries ran at almost 5700 sig/s. The overall slowdown is not due to our query architecture or network latency but more likely because we could not optimize the parameters without direct access to the HSM. Indeed, we asked them to increase the number of threads for the algebraic queries, which more than doubled the signature speed.
Using multiple HSMs as query oracles shows that our attack is mostly platform-agnostic since it only makes black-box use of the HSM functionality for RSA signatures. Even currently used HSMs can be in a configuration state that allows us to run our attack. Furthermore, the Luna S750 HSM test shows that it is possible for an adversary to run these queries over the internet. As shown in Figure 4, we made a total of 4,067,419,331 queries.
4.6 Trip report: Individual signature forgery
The three phases of the individual forgery (descent, linear algebra, and the eth root) all differed in some substantial ways from usual NFS cases. To start off, we generated a random target (t \leftarrow$ Z^*_N).
Descent (i). The first step is descent, which, again, can be broken down into several substeps. The algorithmic ideas are described in Section 3; here we present the progression of the actual 1024-bit computation. Recall from Section 3 that we had four strategies to tackle descent, each specialized to slightly different phases. These were initialization (which closely followed the ideas of [28], and which we can more descriptively call las-plus-ECM); the CADO-NFS GNFS-DL program las descent; composite sieving; and embedded-polynomial sieving. These last two were the main new ideas we applied to descent, as they have previously not been used in the middle and bottom stages. In part, this is because for our chosen smoothness bounds, the difficulty of descent (compared to one-sided precomputation) is unique to the (\tilde{\sqrt{\mathrm{NFS}}}) algorithm.
First, to initialize descent, we greatly benefitted from launching hundreds of jobs, then picking the smoothest (u/v) split we saw. We spent a day running 500 initializations in parallel, each using a different mask and sieving to find (t \cdot \mathrm{mask}^e = u/v \bmod N). We used smoothness parameters lpb0 = lpb1 = 135, very large cofactor bounds mfb0 = mfb1 = 365, sieving area (A = 2^{31}), and factor base bounds lim0 = lim1 = (2^{31}). To perform the large ECM cofactorizations we followed the strategies and scripts of [28]; we printed surviving relations from las to files, on which we ran multithreaded ECM jobs.
After this day of sieving and ECM, we had a clear best initial split: (u/v) factored over primes bounded by 80 bits. In fact, we were only aiming for 135 bits (in the lpb parameters above), but took this lucky split. Because of this, we skipped a step we had needed for the 768-bit computation: the “large q” step we described in Section 3 which performs the jump from 135 to 95 bits. Thus initialization brought us all the way down to 80 bits.
Descent (ii)–(iii). The remainder of descent brought us down from 80 bits (on the rational side) to our final smoothness bounds of (B_{\mathrm{rat}} = 2^{36}) and (B_{\mathrm{ext}} = 2^{35}). In brief, we used the following strategies:
- For most primes from 80 bits to 40 bits: las descent
- For straggler primes from 80 bits to 40 bits: composite sieving or embedded-polynomial sieving, with a few hundred special-q or special-q values
- For primes from 40 bits to the end: composite sieving, with thousands of special-q or special-q values
A strategy which did not work well was simply running las descent with larger and larger A. Even on the machines that could allocate up to (A = 2^{39}), sieving was outperformed by the more parallelizable composite sieving with smaller A.
The parameters we used for the above three strategies were closely tuned to the machines we used; our goal was to use most of the available RAM while maintaining reasonable multi-threading. The memory margin parameter was also critical to ensure we didn’t over-allocate memory. For las descent we used (A = 2^{37}), (\lim = 2^{31}), bkthresh1 = (2^{29}), and bkmult = 1.477. For composite sieving at the bottom, we decreased all these parameters. We used (A = 2^{34}), lim0 = (2^{28}), lim1 = (2^{30}), bkthresh1 = (2^{27}), and bkmult = 1.365.
A final reflection on this 1024-bit descent computation is that although composite sieving worked effectively in practice, it was difficult to predict its success theoretically. Recall that here, when we have a target q to cancel in descent, we sieve over special-q values of q, (q \cdot 11), (q \cdot 13), etc. with small prime (or composite) multipliers. These sieves are not independent, but in practice, certain special-q values seemed to have an easier time in las. Even if one of the first approximately 150 primes showed up in a smooth relation, that relation was sometimes only found when a larger prime multiplier was used. All primes other than the last eight were found (meaning: we found a smooth relation involving that prime) with up to 10,000 small multipliers; most were found early on. For those last eight—which were not remarkable in other ways—we went up to 20,000 small multipliers. Although this tail end of descent was small, it did affect the wall time of descent.
Our ultimate descent tree included 5,633 rational primes bounded by (B_{\mathrm{rat}}), 8,592 prime ideals bounded by (B_{\mathrm{alg}}), and 1,860 prime ideals bounded by (B_{\mathrm{ext}}). These were spread over 874 relations.
Linear algebra. We give a further breakdown of the linear algebra step in Figure 5. Our first observation is that the CADO-NFS bwc program—after some debugging and fixes—handled our linear algebra system quite well. In part, this is because linear algebra in (\tilde{\sqrt{\mathrm{NFS}}}) is an easier task than in GNFS-DL. The modulus (e = 65537) is very small, the factor base of primes (columns) includes only one side, and thanks to extension sieving, we could use an eminently reasonable smoothness bound (B_{\mathrm{alg}} = 2^{31}). On the other hand, SageMath and Python stopped working well (if at all) on matrices of this size. In the 1024-bit run, we omitted several assertions and multiplications since they ran out of memory on our 1TB machine. It seems difficult, in general, to handle a large, sparse, modular matrix.
Despite this, we modified the critical code paths to work. The most important parameter was bwc_mn = 10, which determines how many solutions are found by bwc and how some steps are parallelized. We performed two linear algebra runs: one with additional parallelism via MPI [56], and one without it. The latter instance ran for weeks without interruption and found a correct solution. This is the linear algebra step reported in our timings. The MPI run required substantial manual intervention, but not counting the idle time during breaks, it found a correct answer in shorter wall time but slightly larger CPU time. Concretely, the MPI run used roughly twice the number of machines compared to the non-MPI instance, and had the following timings for the two most compute-intensive steps of linear algebra:
- bwc krylov: 180.73 h (wall time); 80,999.67 h (CPU time); speedup: 1.16x.
- bwc mksol: 179.56 h (wall time); 75,930.18 h (CPU time); speedup: 1.23x.
The other steps had runtimes comparable to the reported non-MPI run. We note that MPI was relatively far away from achieving the maximal possible speedup of 2x in this case. Using additional MPI nodes for every task would further reduce the wall clock time, but eventually see diminishing results due to the communication overhead. Our choice of the number of MPI nodes was due to resource limitations. MPI would likely be more beneficial for larger matrices.
Eth root. Timings of the different root stages are in Figure 5 while the algorithms were explained in Section 3. The parallel Montgomery step could not fully utilize the cores on our machine due to high memory usage, as well as limits to independence in reduction steps. It was also hampered by a few SageMath inefficiencies.
At the start of the root step, our numerator and denominator each had about (2^9) million bits. By the end of the Montgomery step, we had cleared the denominator, and the numerator was bounded, for each coefficient of the number field element, by about (2^{102480}). To get to the very end, our CRT-p-adic step worked well and is included in our timing reports.
We ran a few other tests, not included in our timing. A purely p-adic instance (after Montgomery) with (p = 70111) found a solution after lifting up to (70111^{256}) in about 81 hours of wall time. We ran another Montgomery step with different parameters, and reached an improved numerator bound of (10^{468}). We certainly believe there is room for optimization in the root step, but were satisfied with the reported results.
5 What we talk about when we talk about RSA
The JNT algorithm has been public for two decades, and the research community has been investigating the relationship between RSA and factoring for several more decades. In this section we discuss some consequences (and non-consequences) of (\tilde{\sqrt{\mathrm{NFS}}}) becoming feasible to run on large moduli.
Why doesn’t (\tilde{\sqrt{\mathrm{NFS}}}) break RSA? To start off, we will continue to emphasize the difference between RSA and delayed-target RSA. Although JNT’s paper is titled “When eth roots become easier than factoring” it is more accurately “When eth roots with an oracle become easier than eth roots without an oracle.” In the standard RSA problem, the goal of an attacker, given a public key (N, e) and target t, is to output (t^{1/e} \bmod N). Without the oracle, (\tilde{\sqrt{\mathrm{NFS}}}) cannot be run as is. Even if one could run descent and the rest of the individual forgery phase, it is not obvious how the forgery could be assembled without all the eth roots (coming from oracle queries) on the forgery components in Eq. (1).
One idea is to collect signatures in the wild: for a particular (N, e), could you observe billions of public signatures to replace the oracle queries? In practice, valid signed messages are hashes that are padded to the length of the modulus (albeit often with fixed-pattern padding as in PKCS#1v1.5). These messages will never correspond to the exact rational primes, and rarely to the smooth relations found during precomputation. One might hope to occasionally observe signed values that are smooth with respect to the factor base, and after observing enough such signed values, solve for the signatures on each of the factor base primes. However, this is likely to be much more expensive than the GNFS without additional cleverness. Such a strategy is asking for many full-sized integers to be smooth, which occurs with much smaller probability than the integers sieved over during the GNFS.
Why doesn’t (\tilde{\sqrt{\mathrm{NFS}}}) break factoring? Another observation is the following. The (\tilde{\sqrt{\mathrm{NFS}}}) algorithm allows someone to, after the precomputation, forge any RSA signature of choice. From the outside—and completely ignoring timing—it seems that the attacker must know the secret key, or at least as much information as the secret key.
We do not have any formal arguments, but have two intuitive reasons that (\tilde{\sqrt{\mathrm{NFS}}}) is different from factoring. First, recall that factoring N is equivalent to computing square roots mod N. Even if one had the power to compute arbitrary odd roots mod N; say, (t_1^{1/3}), (t_2^{1/5}), (t_3^{1/7}), and so on for any targets of choice; it is not clear how they could be efficiently combined into an inverse even power.
A second observation we found interesting is that technically, (\tilde{\sqrt{\mathrm{NFS}}}) is only faster than GNFS if the goal is to forge a reasonable number of target signatures. If the goal is to forge subexponentially many signatures for a particular N, for example, it is likely more efficient to run GNFS and learn d than to run the descent, linear algebra, and root steps subexponentially many times. Although we cannot imagine a scenario where this would make sense, there is certainly a range of possible attack strategies.
5.1 One-more and delayed-target RSA
Finally, we explain why we have called our setting the delayed-target RSA problem (DT-RSA), rather than the more popular one-more RSA problem (OM-RSA) [8]. Indeed, it was claimed in [42, §4] that the (\tilde{\sqrt{\mathrm{NFS}}}) breaks a version of the OM-RSA problem. With a cursory analysis, however, this is misleading.
It is true that (\tilde{\sqrt{\mathrm{NFS}}}) solves the DT-RSA problem, and it is true that a polynomial-time attack for DT-RSA yields one for OM-RSA. It is even true that a subexponential-time attack for DT-RSA yields one for OM-RSA. However, when the obvious modifications are applied to (\tilde{\sqrt{\mathrm{NFS}}}) to work in the OM-RSA game, its subexponential parameters are so much worse that they are outperformed by GNFS factoring.
The high-level takeaway of this distinction is that in spirit, (\tilde{\sqrt{\mathrm{NFS}}}) does not attack a one-more problem; it attacks a CCA1-style problem. It is certainly possible that there is a (\tilde{\sqrt{\mathrm{NFS}}})-inspired algorithm that more efficiently attacks OM-RSA, and we do not claim otherwise. The fact that blind signatures tend to rely on the OM-RSA security assumption is interesting. The (\tilde{\sqrt{\mathrm{NFS}}}) algorithm seems to be an eminently relevant and practical attack, but its speedup over factoring is not captured (without modification) by the OM-RSA problem.
An attack on DT-RSA yields an attack on OM-RSA. Formal game-based definitions for the DT-RSA and OM-RSA problems are given in Figure 6. The games are parameterized by key-generation algorithm Kg which for RSA is required to output (N, e, d) such that N = pq is a biprime and (ed = 1 \bmod \phi(N)). They are also parameterized by a number (\ell) of targets, and for OM-RSA, a number r of queries. In the below proposition we prove that an attack on DT-RSA yields one on OM-RSA.
The proposition is not meant to be exhaustively formalized; it is intended to be simple and to illustrate a gap while abstracting the details of (\tilde{\sqrt{\mathrm{NFS}}}). This may lead to better future ideas.
Proposition 1. Let Kg be an RSA key-generation algorithm and let (\ell, r \in \mathbb{N}). Given an adversary B in game (G^{\mathrm{dt-rsa}}{\mathrm{Kg},\ell}) making r queries we can build an adversary A in game (G^{\mathrm{om-rsa}}{\mathrm{Kg},\ell,r}) with
[
\mathrm{Adv}^{\mathrm{dt-rsa}}{\mathrm{Kg},\ell}(B)
\leq
\mathrm{Adv}^{\mathrm{om-rsa}}{\mathrm{Kg},\ell,r}(A).
]
Adversary A also makes r queries. Its running time is about that of B.
Proof. The proof is straightforward. Given adversary B for DT-RSA, we construct adversary A for OM-RSA. A receives as input ((N, e, t_1, \ldots, t_\ell)). It runs B on inputs N and e; the queries to B’s supposed GetQuery oracle are answered by A forwarding the same inputs to and outputs from its own oracle.
When B finishes precomputation and makes its GetTgt query, A returns targets ((t_1, \ldots, t_\ell)). When B finally submits its (\ell) forgeries, A returns the same ones. These are all correct eth roots as long as B was correct. If B made r queries to GetQuery, then so did A. The running time of A is about the same as B. (\square)
Why doesn’t (\tilde{\sqrt{\mathrm{NFS}}}) break OM-RSA? We can use the above proposition to build an attack A on OM-RSA given the attack (\tilde{\sqrt{\mathrm{NFS}}}) on DT-RSA. The main difficulty will be the parameter r, the number of queries.
Throughout this whole paper, we have implicitly assumed the number of targets (\ell) is “reasonable.” (In our computations, (\ell = 1).) When (\tilde{\sqrt{\mathrm{NFS}}}) is run in the (\ell)-target version of DT-RSA, the obvious approach is to run precomputation once (taking time (T_{\mathrm{precomp}})), and the individual forgery (\ell) times (taking time (\ell \cdot T_{\mathrm{indiv}})). Plugging in the relevant costs analyzed by JNT, we have:
[
r = L_N\left[1/3, (8/9)^{1/3} + o(1)\right]
]
[
T_{\mathrm{precomp}} = L_N[1/3, 1.577 + o(1)]
]
[
T_{\mathrm{indiv}} = L_N\left[1/3, (32/9)^{1/3} + o(1)\right].
]
For the minimum nontrivial value of (\ell) ((\ell = r + 1)) in the OM-RSA game, we have that adversary A attacking OM-RSA takes time:
[
T_{\mathrm{precomp}} + \ell \cdot T_{\mathrm{indiv}}
\approx L_N[1/3, 2.488 + o(1)]
]
because (L_N(\alpha, c_1)\cdot L_N(\alpha, c_2)=L_N(\alpha, c_1+c_2)). Recall that the GNFS takes time (L_N(1/3, c + o(1))) with (c = (64/9)^{1/3} \approx 1.923).
Although it could be argued that (\tilde{\sqrt{\mathrm{NFS}}}) thus yields a subexponential-time attack on OM-RSA, that attack is much worse than factoring N in all cases where OM-RSA is nontrivial. The mismatch comes from the query accounting: (\tilde{\sqrt{\mathrm{NFS}}}) requires many queries to be able to forge any single target, while OM-RSA expects (\ell) forgeries from (\ell - 1) queries. One might hope that (\tilde{\sqrt{\mathrm{NFS}}}) could be run with fewer queries, or be modified to take advantage of the OM-RSA adaptive oracle, but these changes would require different algorithmic ideas and parameter analyses (including memory, if parallelism is used on the (\ell) individual forgeries).
5.2 Another type of oracle from JNT’07
As a second revisiting of past work, (\tilde{\sqrt{\mathrm{NFS}}}) does not achieve SNFS complexity when the oracle is restricted to AMR (affine modular root) queries. This difficulty is pointed out in [42, §4.3], and is due to the fact that an AMR oracle imposes a strong constraint on the choice of the algebraic setting of GNFS. The best complexity in that case is given as (L_N(1/3, 6^{1/3} + o(1))) by [42, §4.3], and the (nontrivial) details to achieve that are only given in the ePrint version [43].
Since (6^{1/3} \approx 1.817), we expect that variant to be much less practical than the one discussed in this work, and perhaps even barely competitive with factoring.
The AMR setting best matches the PKCS#1v1.5 signature padding that was until recently ubiquitous in practice for RSA signatures. Unfortunately this asymptotic running time precluded implementing the attack with a more realistic padding oracle given our limited computational resources.
6 Scaling the attack to larger moduli
We can attempt to apply the asymptotic complexity formulas to extrapolate bit securities for larger key sizes than we were able to run, but there are a few difficulties. First, ignoring the hidden constants and o(1) terms in the L notation is a bodacious simplification [46] that can hide real differences in practice; second, we have no guarantees that we have reached asymptotic behavior at 1024 bits; and third, we have no guarantee that our parameter choices for the computations we ran were at all optimal. In fact, during the course of this work, we encountered numerous bugs and infelicities across the software we used—including CADO-NFS, SageMath, and our own code. It is likely that additional issues remain undetected, which could potentially lead to an overestimation of the costs.
With those caveats in place, we extrapolated bit securities of this algorithm for larger key sizes using two methods. In the first method, we computed the multiplicative constant scale factor to recover the generally accepted bit strength for RSA key sizes using the asymptotic running time of the general number field sieve for factoring (blissfully ignoring the o(1) terms), and then applied this constant to the (L_N(1/3, 1.577)) overall running time of our algorithm. Second, we estimated the multiplicative factor for the cycle count of our 1024-bit computation on a 2.3GHz CPU processor and scaled it up to larger bit sizes using the asymptotic complexity formula. Surprisingly, these values were in fairly close alignment, although extrapolation from the smaller computations we did was less effective in predicting our 1024-bit running times.
To estimate the number of queries required, we applied the asymptotic formula for queries to the empirical value for 1024, ignoring o(1) terms and constants. As we discussed in Section 4.3, there is a fair amount of leeway in algorithmic parameter selection, and an attacker could make various tradeoffs to somewhat reduce query complexity, online computation cost, or precomp