Beating the odds on the Diffie-Hellman decision problem

Posted on 30 July 2019 by John

There are a couple variations on the Diffie-Hellman problem in cryptography: the computation problem (CDH) and the decision problem (DDH). This post will explain both and give an example of where the former is hard and the latter easy.

The Diffie-Hellman problems

The Diffie-Hellman problems are formulated for an Abelian group. The main group we have in mind is the multiplicative group of non-zero integers modulo a large prime p. But we start out more generally because the Diffie-Hellman problems are harder over some groups than others as we will see below.

Let g be the generator of a group. When we think of the group operation written as multiplication, this means that every element of the group is a power of g.

Computational Diffie-Hellman problem

Let a and b be two integers. Given g^a and g^b, the CDH problem is to compute g^ab. Note that the exponent is ab and not a+b. The latter would be easy: just multiply g^a and g^b.

To compute g^ab we apparently need to know either a or b, which we are not given. Solving for either exponent is the discrete logarithm problem, which is impractical for some groups.

It’s conceivable that there is a way to solve CDH without solving the discrete logarithm problem, but at this time the most efficient attacks on CDH compute discrete logs.

Diffie-Hellman key exchange

Diffie-Hellman key exchange depends on the assumption that CDH is a hard problem.

Suppose Alice and Bob both agree on a generator g, and select private keys a and b respectively. Alice sends Bob g^a and Bob sends Alice g^b. This doesn’t reveal their private keys because the discrete logarithm problem is (believed to be) hard. Now both compute g^ab and use it as their shared key. Alice computes g^ab by raising g^b to the power a, and Bob computes g^ab by raising g^a to the power b.

If someone intercepted the exchange between Alice and Bob, they would possess g^a and g^b and would want to compute g^ab. This the the CDH.

When working over the integers modulo a large prime p, CDH is hard if p-1 has a large factor, such as when p – 1 = 2q for a prime q. If p-1 has only small factors, i.e. if p-1 is “smooth”, then the discrete logarithm problem is tractable via the Silver-Pohlig-Hellman algorithm [1]. But if p is large and p-1 has a large prime factor, CDH is hard as far as we currently know.

Decision Diffie-Hellman problem

The DDH problem also starts with knowledge of g^a and g^b. But instead of asking to compute g^ab it asks whether one can distinguish g^ab from g^c for a randomly chosen c with probability better than 0.5.

Clearly DDH is weaker than CDH because if we can solve CDH we know the answer to the DDH question with certainty. Still, the DDH problem is believed to be hard over some groups, such as certain elliptic curves, but not over other groups, as illustrated below.

DDH for multiplicative group mod p

Suppose we’re working in the multiplicative group of non-zero integers modulo a prime p. Using Legendre symbols, which are practical to compute even for very large numbers, you can determine which whether g^a is a square mod p, which happens if and only if a is even. So by computing the Legendre symbol for g^a mod p, we know the parity of a. Similarly we can find the parity of b, and so we know the parity of ab. If ab is odd but g^c is a square mod p, we know that the answer to the DDH problem is no. Similarly if ab is even but g^c is not a square mod p, we also know that the answer to the DDH problem is no.

Since half the positive integers mod p are squares, we can answer the DDH problem with certainty half the time by the parity argument above. If our parity argument is inconclusive we have to guess. So overall we can answer he DDH problem correctly with probability 0.75.

[1] As is common when you have a lot of names attached to a theorem, there were multiple discoveries. Silver discovered the algorithm first, but Pohlig and Hellman discovered it independently and published first.

Fame, difficulty, and usefulness

Posted on 24 July 2019 by John

Pierre Fermat is best known for two theorems, dubbed his “last” theorem and his “little” theorem. His last theorem is famous, difficult to prove, and useless. His little theorem is relatively arcane, easy to prove, and extremely useful.

There’s little relation between technical difficulty and usefulness.

Fermat’s last theorem

Fermat’s last theorem says there are no positive integer solutions to

aⁿ + bⁿ = cⁿ

for integers n > 2. This theorem was proven three and a half centuries after Fermat proposed it. The theorem is well known, even in pop culture. For example, Captain Picard tries to prove it in Star Trek: Next Generation. Fermat’s last theorem was famous for being an open problem that was easy to state. Now that it has been proven, perhaps it will fade from popular consciousness.

The mathematics developed over the years in attempts to prove Fermat’s last theorem has been very useful, but the theorem itself is of no practical importance that I’m aware of.

Fermat’s little theorem

Fermat’s little theorem says that if p is a prime and a is any integer not divisible by p, then a^{p − 1} − 1 is divisible by p. This theorem is unknown in pop culture but familiar in math circles. It’s proved near the beginning of any introductory number theory class.

The theorem, and its generalization by Euler, comes up constantly in applications of number theory. See three examples here.

Twisted elliptic curves

Posted on 23 July 2019 by John

This morning I was sitting at a little bakery thinking about what to do before I check out of my hotel. I saw that the name of the bakery was Twist Bakery & Cafe, and that made me think of writing about twisted elliptic curves when I got back to my room.

Twist of an elliptic curve

Suppose p is an odd prime and let E be the elliptic curve with equation

y² = x³ + ax² + bx + c

where all the variables come from the integers mod p. Now pick some d that is not a square mod p, i.e. there is no x such that x² – d is divisible by p. Then the curve E‘ with equation

dy² = x³ + ax² + bx + c

is a twist of E. More explicit it’s a quadratic twist of E. This is usually what people have in mind when they just say twist with no modifier, but there are other kinds of twists.

Let’s get concrete. We’ll work with integers mod 7. Then the squares are 1, 4, and 2. (If that last one looks strange, note that 3*3 = 2 mod 7.) And so the non-squares are 3, 5, and 6. Suppose our curve E is

y² = x³ + 3x + 2.

Then we can form a twist of E by multiplying the left side by 3, 5, or 6. Let’s use 3. So E‘ given by

3y² = x³ + 3x + 2

is a twist of E.

Twisted curve

There’s something potentially misleading in the term “twisted curve”. The curve E‘ is a twist of E, so E‘ is twisted relative to E, and vice versa. You can’t say E‘ is twisted and E is not.

But you might object that E‘ has the form

dy² = x³ + ax² + bx + c

where d is a non-square, and E does not, so E‘ is the one that’s twisted. But a change of variables can leave the curve the same while changing the form of the equation. Every curve has an equation in which the coefficient of y² is 1 and a form in which the coefficient is a non-square. Specifically,

dy² = x³ + ax² + bx + c

and

y² = x³ + dax² + d²bx + d³c

specify the same curve.

The two curves E and E‘ are not isomorphic over the field of integers mod p, but they are isomorphic over the quadratic extension of the integers mod p. That is, if you extend the field by adding the square root of d, the two curves are isomorphic over that field.

Partitioning x coordinates

For a given x, f(x) = x³ + ax² + bx + c is either a square or a non-square. If it is a square, then

y² = f(x)

has a solution and

dy² = f(x)

does not. Conversely, if f(x) is not a square, then

dy² = f(x)

has a solution and

y² = f(x)

does not. So a given x coordinate either corresponds to a point on E or a point on E‘, but not both.

Application to cryptography

In elliptic curve cryptography, it’s good if not only is the curve you’re using secure, so are its twists. Suppose you intend to work over a curve E, and someone sends you an x coordinate that’s not on E. If you don’t check whether x corresponds to a point on E, you are effectively working on a twist E‘ rather than E. That can be a security hole if E‘ is not as secure a curve as E, i.e. if the discrete logarithm problem on E‘ can be solved more efficiently than the same problem on E.

Distribution of quadratic residues

Posted on 12 July 2019 by John

Let p be an odd prime number. If the equation

x² = n mod p

has a solution then n is a square mod p, or in classical terminology, n is a quadratic residue mod p. Half of the numbers between 0 and p are quadratic residues and half are not. The residues are distributed somewhat randomly. For large p, quadratic residues are distributed enough like random numbers to be useful in cryptographic applications. For example, see this post. But the distribution isn’t entirely random as we’ll demonstrate here.

First let’s look at a couple illustrations, using p = 101 and p = 103. We’ll use the following Python code to plot the residues.

    import matplotlib.pyplot as plt

    for p in [101, 103]:
        s = {x**2 % p for x in range(1, p)}
        for q in s:
            plt.vlines(q, 0, 1)
        plt.title("Quadratic residues mod {}".format(p))
        plt.savefig("residues_mod_{}.svg".format(p))
        plt.close()

Here are the plots.

quadratic residues mod 101

quadratic residues mod 103

Symmetry and anti-symmetry

If you look closely at the first image, you’ll notice it appears to be symmetric about the middle. And if you look closely at the second image, you’ll notice it appears to be anti-symmetric, i.e. the left side has bars where the right side has gaps and vice versa.

The following code shows that these observations are correct.

    p = 101
    s = {x**2 % p for x in range(1, p)}
    for q in s:
        assert(p - q in s)
    
    p = 103
    s = {x**2 % p for x in range(1, p)}
    for q in s:
        assert(p - q not in s)

Both of these observations hold more generally. If p = 1 mod 4, the residues are symmetric: n is a quadratic residue if and only if p – n is a quadratic residue. And if p = 3 mod 4, the residues are anti-symmetric: n is a quadratic residue if and only if p – n is not a quadratic residue. Both of these statements are consequences of the quadratic reciprocity theorem.

Quadratic excess

If you look again at the plot of residues for p = 103, you’ll notice that not only is it anti-symmetric, it is also darker on the left side. We can verify with the following code

   print(len( [q for q in s if q < 51] )) 
   print(len( [q for q in s if q > 51] ))

that there are 28 residues on the left and and 23 on the right. This is a general pattern called quadratic excess. The quadratic excess of an odd prime p, denoted E(p) is the number of quadratic residues to the left of the middle minus the number to the right of the middle. If p = 1 mod 4, E(p) is zero by symmetry, but if p = 3 mod 4, E(p) is always positive.

Here is a plot of quadratic excess for primes less than 3000 and congruent to 3 mod 4.

plot of quadratic excess

Statistical properties

In this post I’ve said a couple things that are in tension. At the top I said that quadratic residues act sufficiently like random bits to be useful in cryptography. Then later on I showed that they have major departures from random behavior. That would be a problem if you looked at the entire sequence of residues, but if p has thousands of digits, and you’re looking at a small sample of the reciprocity values, that’s a different situation.

Suppose you pick two large enough primes, p and q, and keep them secret but tell me their product N = pq. Then you pick a random value n and ask me whether it’s a quadratic residue mod N, the best I can do is say “Maybe, with probability 1/2.” There are crypto systems that depend on this being a computationally hard problem to solve, unless the factors p and q are known.

However, you need to be a little careful. If n is smaller than √N, I could test whether n is a square, working in the integers, not the integers mod N. If it is a square integer, then it’s a square mod N. If it’s not a square integer, it might be a square mod N, but I wouldn’t know.

There’s a related issue I may explore in the future. Suppose instead of just asking whether a particular number is a quadratic residue, you take a range of numbers and return 0 or 1 depending on whether each number is a residue. Can I distinguish the sequence from a random sequence using statistical tests? If I took a sequence of length greater than N/2 then I could tell by looking at symmetry or antisymmetry. But what if the range is small relative to N? Maybe N is a thousand-digit numnber but you only give me a sequence of a million bits. Would the distribution be distinguishable from uniform random bits? Is there any sort of autocorrelation that might betray the sequence as not random?

Software to factor integers

Posted on 25 June 2019 by John

In my previous post, I showed how changing one bit of a semiprime (i.e. the product of two primes) creates an integer that can be factored much faster. I started writing that post using Python with SymPy, but moved to Mathematica because factoring took too long.

SymPy vs Mathematica

When I’m working in Python, SymPy lets me stay in Python. I’ll often use SymPy for a task that Mathematica could do better just so I can stay in one environment. But sometimes efficiency is a problem.

SymPy is written in pure Python, for better and for worse. When it comes to factoring large integers, it’s for worse. I tried factoring a 140-bit integer with SymPy, and killed the process after over an hour. Mathematica factored the same integer in 1/3 of a second.

Mathematica vs PARI/GP

The previous post factors 200-bit semiprimes. The first example, N = pq where

p = 1078376712338123201911958185123
q = 1126171711601272883728179081277

took 99.94 seconds to factor using Mathematica. A random sample of 13 products of 100-bit primes and they took an average of 99.1 seconds to factor.

Using PARI/GP, factoring the value of N above took 11.4 seconds to factor. I then generated a sample of 10 products of 100-bit primes and on average they took 10.4 seconds to factor using PARI/GP.

So in these examples, Mathematica is several orders of magnitude faster than SymPy, and PARI/GP is one order of magnitude faster than Mathematica.

It could be that the PARI/GP algorithms are relatively better at factoring semiprimes. To compare the efficiency of PARI/GP and Mathematica on non-semiprimes, I repeated the exercise in the previous post, flipping each bit of N one at a time and factoring.

This took 240.3 seconds with PARI/GP. The same code in Mathematica took 994.5 seconds. So in this example, PARI/GP is about 4 times faster where as for semiprimes it was 10 times faster.

Python and PARI

There is a Python interface to PARI called cypari2. It should offer the convenience of working in Python with the efficiency of PARI. Unfortunately, the installation failed on my computer. I think SageMath interfaces Python to PARI but I haven’t tried it.

Making public keys factorable with Rowhammer

Posted on 21 June 2019 by John

The security of RSA encryption depends on the fact that the product of two large primes is difficult to factor. So if p and q are large primes, say 2048 bits each, then you can publish n = pq with little fear that someone can factor n to recover p and q.

But if you can change n by a tiny amount, you may make it much easier to factor. The Rowhammer attack does this by causing DRAM memory to flip bits. Note that we’re not talking about breaking someone’s encryption in the usual sense. We’re talking about secretly changing their encryption to a system we can break.

To illustrate on a small scale what a difference changing one bit can make, let p = 251 and q = 643. Then n = pq = 161393. If we flip the last bit of n we get m = 161392. Although n is hard to factor by hand because it has no small factors, m is easy to factor, and in fact

161392 = 2⁴ × 7 × 11 × 131.

For a larger example, I generated two 100-bit random primes in Mathematica

p = 1078376712338123201911958185123
q = 1126171711601272883728179081277

and was able to have it factor n = pq in about 100 seconds. But Mathematica was able to factor n-1 in a third of a second.

So far we have looked at flipping the least significant bit. But Rowhammer isn’t that precise. It might flip some other bit.

If you flip any bit of a product of two large primes, you’re likely to get an easier factoring problem, though the difficulty depends on the number you start with and which bit you flip. To illustrate this, I flipped each of the bits one at a time and measured how long it took to factor the result.

The median time to factor n with one bit flipped was 0.4 seconds. Here’s a plot of the factoring times as a function of which bit was flipped.

factoring times for corrupted product of primes

The plot shows about 80% of the data. Twenty percent of the time the value was above 11 seconds, and the maximum value was 74 seconds. So in every case flipping one bit made the factorization easier, usually quite a lot easier, but only a little easier in the worst case.

To verify that the results above were typical, I did a second experiment. This time I generated a sequence of pairs of random 100-bit primes. I factored their product, then factored the product with a randomly chosen bit flipped. Here are the factoring times in seconds.

    |----------+---------|
    | Original | Flipped |
    |----------+---------|
    |  117.563 |   3.828 |
    |  108.672 |   4.875 |
    |   99.641 |   0.422 |
    |  103.031 |   0.000 |
    |   99.188 |   0.000 |
    |  102.453 |   0.234 |
    |   79.594 |   0.094 |
    |   91.031 |   0.875 |
    |   64.313 |   0.000 |
    |   95.719 |   6.500 |
    |  131.125 |   0.375 |
    |   97.219 |   0.000 |
    |   98.828 |   0.203 |
    |----------+---------|

By the way, we started this post by talking about maliciously flipping a bit. The same thing can happen without a planned attack if memory is hit by a cosmic ray.

Bounds on the nth prime

Posted on 20 June 2019 by John

The nth prime is approximately n log n.

For more precise estimates, there are numerous upper and lower bounds for the nth prime, each tighter over some intervals than others. Here I want to point out upper and lower bounds from a dissertation by Christian Axler on page viii.

First, define

$\begin{align*} f(n, k) &= \log n + \log \log n - 1 + \frac{\log \log n - 2}{\log n} \\ &\qquad - \frac{(\log\log n)^2 - 6\log\log n + k}{2\log^2 n} \end{align*}$

Then for sufficiently large n, the nth prime number, p_n, is bounded above and below by

$n f(n, 11.847) <\, p_n < n f(n, 10.273)$

The lower bound holds for n ≥ 2, and the upper bound holds for n ≥ 8,009,824.

For example, suppose we want to estimate the 10 millionth prime. The exact value is 179,424,673. The bounds we get from the equations above are 179,408,030 and 179,438,842.

The width of the bracket bounding p_n is 0.787 n / log²n. In our example, the difference between the upper and lower bounds is 30,812.

This width, relative to p_n, decreases something like O(1/log³ n). In our example, the width of the bounding interval is 0.017% of p_n.

Feller-Tornier constant

Posted on 18 June 2019 by John

Here’s kind of an unusual question: What is the density of integers that have an even number of prime factors with an exponent greater than 1?

To define the density, you take the proportion up to an integer N then take the limit as N goes to infinity.

It’s not obvious that the limit should exist. But if it does exist, you might guess that it’s 1/2; it’s a question involving whether something is even or odd, and so even and odd occurrences might balance each other out.

But the result, known as the Feller-Tonier constant, is bigger than 1/2. This seems more reasonable when you consider that zero is an even number and a lot of numbers may not have any prime factors with exponent bigger than 1.

Here’s a little Python code to compute the ratio for N = 10⁶.

    from sympy.ntheory import factorint

    def f(n):
        v = factorint(n).values()
        n = len([x for x in v if x > 1])
        return n % 2 == 0

    c = 0
    N = 1000000
    for n in range(N):
        if f(n):
            c += 1
    print (c/N)

This gives a ratio of 0.661344. The limiting value is 0.661370494….

Computing the Feller-Tornier constant

The code above gives an estimate for the Feller-Tornier constant, but how accurate is it? Is there a more efficient way to estimate it that doesn’t require factoring a lot of numbers?

The constant is given analytically by

$c = \frac{1}{2} + \frac{1}{2}\prod_{n=1}^\infty \left(1 - \frac{2}{p_n^2} \right )$

where p_n is the nth prime.

Here’s some Python code to estimate the Feller-Tornier constant using the expression above.

    N = 1000000
    product = 1
    for p in primerange(2, N):
        product *= 1 - 2*p**-2
    print(0.5*(1+product))

Note that we used primerange rather than prime(n). If you need to generate a range of prime numbers, the former is much more efficient.

The code gives us an estimate of 0.66131707, which is good to seven decimal places. Furthermore, we know this is an upper bound on the exact value because we’ve left out terms in the infinite product that are less than 1.

One could calculate a lower bound as well and thus know an interval that must contain the true value. The lower bound is left as an exercise to the reader.

Strong primes

Posted on 20 May 2019 by John

There are a couple different definitions of a strong prime. In number theory, a strong prime is one that is closer to the next prime than to the previous prime. For example, 11 is a strong prime because it is closer to 13 than to 7.

In cryptography, a strong primes are roughly speaking primes whose products are likely to be hard to factor. More specifically, though still not too specific, p is a strong prime if

p is large
p – 1 has a large prime factor q
q – 1 has a large prime factor r
p + 1 has a large prime factor s

The meaning of “large” is not precise, and varies over time. In (1), large means large enough that it is suitable for use in cryptography, such as in RSA encryption. This standard increases over time due to increases in computational power and improvements in factoring technology. The meaning of “large” in (2), (3), and (4) is not precise, but makes sense in relative terms. For example in (2), the smaller the ratio (p – 1)/q the better.

Relation between the definitions

The Wikipedia article on strong primes makes the following claim without any details:

A computationally large safe prime is likely to be a cryptographically strong prime.

I don’t know whether this has been proven, or even if it’s true, but I’d like to explore it empirically. (Update: see the section on safe primes below. I misread “safe” above as “strong.” Just as well: it lead to an interesting investigation.)

We’ll need some way to quantify whether a prime is strong in the cryptographic sense. This has probably been done before, but for my purposes I’ll use the sum of the logarithms of q, r, and s. We should look at these relative to the size of p, but all the p‘s I generate will be roughly the same size.

Python code

I’ll generate 100-bit primes just so my script will run quickly. These primes are too small for use in practice, but hopefully the results here will be representative of larger primes.

    from sympy import nextprime, prevprime, factorint, randprime
    import numpy as np
    
    # largest prime factor
    def lpf(n):
        return max(factorint(n).keys())
    
    def log2(n):
        np.log2(float(n))
    
    num_samples = 100
    data = np.zeros((num_samples, 5))
    
    bitsize = 100
    
    for i in range(num_samples):
        p = randprime(2**bitsize, 2**(bitsize+1))
        data[i,0] = 2*p > nextprime(p) + prevprime(p)
        q = lpf(p-1)
        r = lpf(q-1)
        s = lpf(p+1)
        data[i,1] = log2(q)
        data[i,2] = log2(r)
        data[i,3] = log2(s)
        data[i,4] = log2(q*r*s)

The columns of our matrix correspond to whether the prime is strong in the number theory sense, the number of bits in q, r, and s, and the total bits in the three numbers. (Technically the log base 2 rather than the number of bits.)

Results

There were 75 strong primes and 25 non-strong primes. Here were the averages:

    |-----+--------+------------|
    |     | strong | not strong |
    |-----+--------+------------|
    | q   |   63.6 |       58.8 |
    | r   |   41.2 |       37.0 |
    | s   |   66.3 |       64.3 |
    | sum |  171.0 |      160.1 |
    |-----+--------+------------|

The numbers are consistently higher for strong primes. However, the differences are small relative to the standard deviations of the values. Here are the standard deviations:

    |-----+--------+------------|
    |     | strong | not strong |
    |-----+--------+------------|
    | q   |   20.7 |       15.6 |
    | r   |   19.8 |       12.3 |
    | s   |   18.7 |       19.9 |
    | sum |   30.8 |       41.9 |
    |-----+--------+------------|

Safe primes

I realized after publishing this post that the Wikipedia quote didn’t say what I thought it did. It said that safe primes are likely to be cryptographically strong primes. I misread that as strong primes. But the investigation above remains valid. It shows weak evidence that strong primes in the number theoretical sense are also strong primes in the cryptographic sense.

Note that safe does not imply strong; it only implies the second criterion in the definition of strong. Also, strong does not imply safe.

To test empirically whether safe primes are likely to be cryptographically strong, I modified my code to generate safe primes and compute the strength as before, the sum of the logs base 2 of q, r, and s. We should expect the strength to be larger since the largest factor of p will always be as large as possible, (p – 1)/2. But there’s no obvious reason why r or s should be large.

For 100-bit safe primes, I got an average strength of 225.4 with standard deviation 22.8, much larger than in my first experiment, and with less variance.

Golden ratio primes

Posted on 12 May 2019 by John

The golden ratio is the larger root of the equation

φ² – φ – 1 = 0.

By analogy, golden ratio primes are prime numbers of the form

p = φ² – φ – 1

where φ is an integer. To put it another way, instead of solving the equation

φ² – φ – 1 = 0

over the real numbers, we’re looking for prime numbers p where the equation can be solved in the integers mod p. [1]

Application

When φ is a large power of 2, these prime numbers are useful in cryptography because their special form makes modular multiplication more efficient. (See the previous post on Ed448.) We could look for such primes with the following Python code.

    from sympy import isprime

    for n in range(1000):
        phi = 2**n
        q = phi**2 - phi - 1
        if isprime(q):
            print(n)

This prints 19 results, including n = 224, corresponding to the golden ratio prime in the previous post. This is the only output where n is a multiple of 32, which was useful in the design of Ed448.

Golden ratio primes in general

Of course you could look for golden ratio primes where φ is not a power of 2. It’s just that powers of 2 are the application where I first encountered them.

A prime number p is a golden ratio prime if there exists an integer φ such that

p = φ² – φ – 1

which, by the quadratic theorem, is equivalent to requiring that m = 4p + 5 is a square. In that case

φ = (1 + √m)/2.

Here’s some code for seeing which primes less than 1000 are golden ratio primes.

    from sympy import primerange

    def issquare(m):
        return int(m**0.5)**2 == m

    for p in primerange(2, 1000):
        m = 4*p + 5
        if issquare(m):
            phi = (int(m**0.5) + 1) // 2
            assert(p == phi**2 - phi - 1)
            print(p)

By the way, there are faster ways to determine whether an integer is a square. See this post for algorithms.

(Update: Aaron Meurer pointed out in the comments that SymPy has an efficient function sympy.ntheory.primetest.is_square for testing whether a number is a square.)

Instead of looping over primes and testing whether it’s possible to solve for φ, we could loop over φ and test whether φ leads to a prime number.

    for phi in range(1000):
        p = phi**2 - phi - 1
        if isprime(p):     
            print(phi, p)

Examples

The smallest golden ratio prime is p = 5, with φ = 3.

Here’s a cute one: the pi prime 314159 is a golden ratio prime, with φ = 561.

The golden ratio prime that started this rabbit trail was the one with φ = 2²²⁴, which Mike Hamburg calls the Goldilocks prime in his design of Ed448.

[1] If p = φ² – φ – 1 for some integer φ, then φ² – φ – 1 = 0 (mod p). But the congruence can have a solution when p is not a golden ratio prime. The following code shows that the smallest example is p = 31 and φ = 13.

from sympy import primerange
from sympy.ntheory.primetest import is_square

for p in primerange(2, 100):
    m = 4*p + 5
    if not is_square(m):
        for x in range(p):
            if (x**2 - x - 1) % p == 0:
                print(p, x)
                exit()

The Diffie-Hellman problems

Computational Diffie-Hellman problem

Diffie-Hellman key exchange

Decision Diffie-Hellman problem

DDH for multiplicative group mod p

Related posts

Fermat’s last theorem

Fermat’s little theorem

Related posts

Twist of an elliptic curve

Twisted curve

Partitioning x coordinates

Application to cryptography

Related posts

Symmetry and anti-symmetry

Quadratic excess

Statistical properties

Related posts

SymPy vs Mathematica

Mathematica vs PARI/GP

Python and PARI

Related posts

Related posts

Computing the Feller-Tornier constant

Relation between the definitions

Python code

Results

Safe primes

Related posts

Application

Golden ratio primes in general

Examples

Related posts