Cyclic groups and cryptography

Elliptic curve cryptography is based upon a cyclic group, in which we consider a "hard" discrete logarithm problem.  As outlined before, the discrete logarithm problem stems from the fact that in a cyclic group with generator G, it is easy to calculate A = Gn where n is a natural number smaller than the order of the group, but it is sometimes a very hard problem to calculate n from A, n being the "discrete logarithm in base G" of the group element A in the cyclic group.  We need a cyclic group of course, because every element in the group is going to be reached by a power of the generator G.  In fact, the denomination of "logarithm" and the use of "exponent" come from the idea that the group has a "multiplicative" law.  But if our group has an "additive" law, then our "exponent" simply becomes a "multiplication", and our "logarithm" becomes a "division".  So then we have the "hard division problem".  But not all division problems are hard !  We don't want our group to be the "additive group" of a field.  Indeed, in a field, the division is easily defined, and would solve the "discrete division problem" trivially. The historical reason why we thought of a "multiplicative" group, was that it is algebraically known that we cannot have "triple fields", that is, a set with an additive, multiplicative, and exponential law so that addition and multiplication form a field, and that multiplication and exponentiation form another field.  As such, it was considered "safe" to consider the multiplicative group of an existing field.  As it turns out that the multiplicative group of a prime modular field is cyclic, this is a good candidate to have a discrete logarithm problem in it, and it is used a lot.  However, the discrete logarithm problem in those groups has been studied quite a lot, and a lot of progress has been made in solving them, so one needs to go to very big groups in order for the discrete logarithm problem to be hard enough to be safe cryptographically.  In fact, the most general attack (general number sieve index method) on integer factorisation as well as discrete logarithm problems in the multiplicative group of a modular field is given by the following formula:

nsecurity = (1.923*(b*ln(2))1/3 ( ln(b*ln(2)) )2/3 / ln(2)

where b is the order of the group and n is the bit security level.  Grossly, we see that the security level goes roughly as b1/3.

This gives rise to the standard list of needed group orders  (logof them, the number of bits in the "key") in order to achieve a given security level:

security level (bits) group order (bits)
80 1024
112 2048
128 3072
192 7680
256 15360

For high levels of security, the group orders (and key lengths) become quite large.

But one can construct other groups, and the family of elliptic groups is such another series of groups.   Algebraically, all cyclic groups of order n are isomorphic, and isomorphic to Zn, +.  So essentially, no matter what cyclic group one is trying to use, one is in fact just using an isomorphism f between Zn, + and A,*, in such a way that f(1) = G, the generator, and f(n) = Gn.  As the group elements in A are to be represented also by integers, h(U) = m, we have essentially that our isomorphism comes down to a complicated numerical function g(n) = h(f(n)) = h(Gn) = m, so that:

  • if we are given k = g(n) and l = g(m), then we can calculate, from k and l alone and the group peculiarity, p, such that p = g(n+m).
  • if we are given k = g(n.m) and m, however, we cannot calculate q = g(n) from k alone (hard discrete logarithm problem), even though if we know r = n.m, and m, there's of course no difficulty in finding n, namely by dividing r by m.

In other words, even though the isomorphism f induces, on the group G, also a ring and even a field structure, we only know the group operation of addition, and we have no idea how to implement the ring division without "going back" and f is considered a trap door function. 

So the actual "cryptographic" element is not the group structure itself, but rather the intractability of the isomorphism between Zn,+ and a particular group.  It is believed that elliptic groups have more intractable isomorphisms than the multiplicative groups of modular fields, even though all of them are isomorphic to Zn, +.

Definition of elliptic groups

Elliptic curves

Elliptic groups are induced by a "discretization" of real elliptic curves.  Elliptic curves are simply curves with the following equation:

y2 = x3 + a x + b

in the same way that a straight line is defined by:

y = a x + b

that a circle is defined by:

y2 = R2 - x2

and that an ellipse is defined by:

y2 = R2 - a x2

in other words, it is just a curve in the plane with a given equation. In fact, there is usually an extra requirement, that the parameters a and b satisfy the condition that:

4 a3 + 27 b2 is NOT equal to zero.

In that equation, there are two parameters, a and b, that determine the exact form of the elliptic curve.  However, this curve has some peculiar properties:

  1. if we take two arbitrary points on the curve, and draw a straight line, then this straight line intersects the curve in a third point if the line is not vertical.
  2. the curve is symmetric in the x-axis: for every point with a positive y, there is also exactly another point with a negative y.

This may seem like a banal property, but in fact, this allows one to associate, to each pair of points on the curve, a third point, which is the y-flipped intersection point of the straight line through the first two.  In other words, if we have two points P and Q, of the curve, then we can associate a third point, R = s(P,Q).  It is pretty obvious that s is a symmetric relationship: s(P,Q) = s(Q,P), because the straight line through P and Q is also the straight line through Q and P.

If we define now the "opposite" of a point R, called -R, the point with flipped y value, we see that we have:

  • the straight line through P and Q hits -R, so R = s(P,Q)
  • the straight line through -R and P hits Q, so -Q = s(P,-R)
  • the straight line through -R and Q hits P, so -P = s(Q,-R)

Note that this is exactly what we would expect from a binary, commutative operation "+" instead of s():

  • P+Q = R
  • P + (-R) = -Q
  • Q + (-R) = -P

The thing we still have to introduce is the neutral element.  It is not a point on the curve, and we introduce a new point "zero" which is the result of s if we have a vertical line, and there's no genuine intersection point with the curve.  So s(R, -R) doesn't exist on the curve, and we say that it is zero, Z.  In the same way, s(R,Z) = R, etc...

As such, we have an internal composition law of points on the curve that is:

  • internal
  • has an artificially introduced neutral element
  • has an inverse element for every element

The remarkable property of elliptic curves, is moreover that associativity hold !

As such, this composition law is a commutative group law !

The algebraic calculation of the (flipped) intersection points  is as follows:

a) if xP is different from xQ: Calculate s = (yP - yQ) / (xp - xQ), and then calculate: xR = s2 - xP - xQ and yR = - yP - s(xR - xP)

b) if xP = xQ and if yP = - yQ then the result is Zero

c) if xP = xQ and if yP = yQ != 0 (so P = Q)  then we calculate s = (3 xP2 + a)/(2 yP) and xR = s2 - 2 xP and yR = - yP - s(xR - xP)

d) if xP = xQ and if yP = yQ = 0 (so P = Q) then the result is the neutral element

e) if one of the elements (Q or P) is the neutral element, the sum is equal to the other element (P or Q)

As the expression of this equation is purely algebraic, and can be done with additions and multiplications only, this equation and the definition  can just as well apply in any field, not only in the field of real numbers.  As the expressions for the initially geometric properties are also purely algebraic, they hold up just as well in any field ; the properties of the composition law thus defined, are then also inducing a group law.

Discrete elliptic groups over a finite field F

Consider a field F, and consider couples of elements in F x F.   Consider those couples that satisfy an elliptic curve equation with chosen parameters a and b in F.  This determines a set of couples, part of F x F.  Add an extra element to this set, which we call "zero", and we call this extended set E.  The composition law as defined above, makes that this extended set is a commutative group: the elliptic group over F, with the chosen parameters a and b.

If F is a finite field, and E is an elliptic group over F with parameters a and b, then one can wonder how many elements this group E will contain.  Hasse's theorem states that if q is the order of the field F, then the elliptic group must have an order n that is not very different from q: the elliptic group will deviate at most 2 sqrt(q) from the order q.  So if q is very large, then the elliptic group is also of order not very far from q.   

Cyclic subgroups

For cryptographic purposes, we don't just need a group ; we need a "complicated isomorphism" between Zp,+ and a cyclic subgroup of the group, with p large.  In a given finite group, one can always construct cyclic subgroups, by simply picking an element of the group and constructing all its multiples / powers (depending whether one considers additive language or multiplicative language).  But the problem is that we need to be sure that our p is very large, and that the isomorphism is intractable, that is, that it doesn't have mathematically known properties that allow one to reverse it easily.

What is interesting, is to know that an elliptic group E over a finite field is a cyclic group, or the product of two cyclic groups (attention: these cyclic groups can have sub-groups themselves of course !).   By picking an element of the group, and by considering all its multiples, one will of course construct a cyclic subgroup of the elliptic group.  However, the danger in doing so is that one has an elliptic group that contains two cyclic subgroups, and that one has picked an element in a small subgroup !

Cyclic subgroup of an elliptic group, and cryptography

Suppose that we have chosen a field F, an elliptic curve (coefficients a and b in F) which defines an elliptic group E(F,a,b), and a generator g which defines a cyclic subgroup G in E(F,a,b).   In order to be secure, we need the cyclic group to be of prime order, so that it doesn't contain small sub-groups.  There is in fact only one solution: that G is E(F,a,b), which happens when the order of E(F,a,b), N, is prime itself.   If the order of E(F,a,b) is prime, then any generator will do, because all elements that are not the neutral element, are of order, the prime order of the group.

A priori, the safety level in bits of a cyclic group of prime order N is (at most) of the order of 1/2 log2 N, because all of these discrete logarithm cases can be attacked by a variant of the Pollard's rho attack, who reduces the number of trials statistically from N to sqrt(N).  It is believed that most elliptic groups have about this safety level, because for most elliptic groups, no more powerful method is known.  This is the main motivation to use elliptic groups, that the security level is thought to be of the level of log(sqrt(N)), and not, like modular field discrete logarithm or factorisation problems, of the level of (log(N))1/3.

So it seems that it is OK to find an elliptic curve with a prime order N over a finite field F.

However, that is not sufficient.  Indeed, reductions of the elliptic group to the multiplicative group of a modular field Fq are known in the following cases:

  • when N = p, the order of the finite field over which the elliptic group was defined
  • when there exists a prime number r, such that (rt - 1) is a multiple of N, then q = rt, so this q should be large enough to give rise to a harder discrete log problem than the one in the elliptic group itself

As such, the construction of a cryptographically reliable discrete logarithm problem with an elliptic group demands some precautions in the choice of the field Fp, and the curve (parameters a and b), determining the (prime) order N of the elliptic group satisfying that N is different from p, and that the only prime powers rt that are multiples of N, plus 1, are large enough.  If one is not very sure about these properties, it is better to use a standard curve than inventing one one-self, but of course, in that case, one needs to trust the body that has provided the standard curve not to possess a back door somehow to that group.

In a way, this is somewhat scary: one has discovered two families of "weak" elliptic curves, where the discrete logarithm problem in the elliptic group can (easily) be transformed in a discrete logarithm problem of a modular multiplicative group, of which it is known that they can be solved with a security level of only (log(N))1/3.  Maybe there are other such families that are still undiscovered, or simply not publicly known !  The only reassuring part is that these families have been discovered 30 years ago, and that 30 years of research hasn't brought in yet another weak family, at least, publicly. 

The idea is that if no other weak families have been found in 30 years of research, that if ever such families exist, they mustn't be too densely present in the set of "good" (until now) candidate curves.  In order to increase the probability to use a non-weak curve, the proposition is to use "randomly" selected curves, not doctored by one or other mathematical criteria, which could exactly be the key to another weak family.  However, in order to be entirely clear, this "randomly" needs to be provable randomly, and not obfuscated doctored.

This last part is a critical aspect in some standardised public curves: the randomness of the "random choice" is not proven.  As such, the possibility exists that these "randomly selected" curves are in fact very doctored curves with some or other not publicly known weakness to it.  On the other hand, several of those standard curves have been published about 20 years ago, and no extra "weak family" has been discovered publicly since that time, so the chances that people cooking those standards had a knowledge of a non-publicly known weak family become smaller and smaller as time goes by.   If one is really paranoid about such possibilities, one can generate one's own random elliptic group, but one needs to take the precautions listed above.

Code: basic elliptic curve operations

We recall the general framework concerning Code examples

For starters, we will not generate ourselves a good elliptic curve with generator, because this involves counting the elements N in the cyclic group obtained, which is quite involved for large groups.  So in the code examples that follow, we will assume that we have a given elliptic group (a given prime, a given curve, and a given generator).

We define two new types:

/* "elliptic" is a structure that contains a number and a boolean that tells us whether
   the number is "zero" or not, because elliptic curves have a non-number neutral element */
typedef struct digitos {bool zero ; number x ; number y;} elliptic;

/* "ellipcurve" is a structure that contains the parameters of an elliptic curve on a modular field */
typedef struct digitossos {number prime ; number a ; number b ; elliptic generator; } ellipcurve;

An elliptic point is a record containing two "numbers", but also a flag that tells us whether the point is the neutral element (not represented by any "numbers" or not. An elliptic curve is specified by the modular prime (we limit ourselves here to fields which are prime fields), the two coefficients a and b of the elliptic curve equation, and, as we will see, the generator of the cyclic subgroup in which we will consider the discrete logarithm problem.

We need an auxiliary function that copies an elliptic point into another one:

void ellipcopy(elliptic* orig, elliptic* copy)
/* copies an elliptic point into another one */
  {
  int i;
  copy->zero = orig->zero;
  for (i = 0 ; i < nofdigits ; i++)
    {
    copy->x.digit[i] = orig->x.digit[i];
    copy->y.digit[i] = orig->y.digit[i];
    }
  }  // end ellipcopy 

The core algorithm which makes the sum of two elliptic points goes as follows:

void ellipsum(ellipcurve* mycurve,elliptic* a, elliptic* b, elliptic* sum)
/* calculates the sum of two elliptic points on a given curve */
  {
  number s;
  number zeero;
  number two;
  number three;
  number xa2;
  number xa2t3;
  number xa2t3pa;
  number yat2;
  number invyat2;
  number xamxs;
  number stxamxs;
  number s2;
  number xat2;
  number yamyb;
  number xamxb;
  number invxamxb;
  number s2mxa;
  number oldmodulo;
  unsigned int i;
  
  // we apply the prime number of the modular field of the curve
  for (i = 0 ; i < nofdigits ; i++)
    {
    oldmodulo.digit[i] = modulo.digit[i];
    modulo.digit[i] = mycurve->prime.digit[i];
    }
  makebig(0,&zeero);
  if (a->zero && b->zero)
    {
    // we have zero plus zero
    sum->zero = true;
    }
  else if (a->zero)
    {
    // a is zero, so the sum is b
    ellipcopy(b,sum);
    }
  else if (b->zero)
    {
    // b is zero, so the sum is a
    ellipcopy(a,sum);
    }
  else if (isbigger(&(a->x),&(b->x)) == 0)
    {
    // both points a and b are on a vertical line
    // there are only two possibilities: a = b != 0 or a = -b
    if (isbigger(&(a->y),&(b->y)) == 0 && isbigger(&(a->y),&zeero) != 0)
      {
      // a is b, so we have 2a
      sum->zero = false;
      timesnum(&(a->x),&(a->x),&xa2);
      makebig(2,&two);
      makebig(3,&three);
      timesnum(&three,&xa2,&xa2t3);
      sumnum(&xa2t3,&(mycurve->a),&xa2t3pa);
      timesnum(&two,&(a->y),&yat2);
      inversenum(&yat2,&invyat2);
      timesnum(&xa2t3pa,&invyat2,&s);
      timesnum(&s,&s,&s2);
      timesnum(&two,&(a->x),&xat2);
      difnum(&s2,&xat2,&(sum->x));
      difnum(&(a->x),&(sum->x),&xamxs);
      timesnum(&xamxs,&s,&stxamxs);
      difnum(&stxamxs,&(a->y),&(sum->y));
      }
    else
      {
      // a is minus b, so we have zero
      sum->zero = true;
      }
    }
  else
    {
    // a and b are different, not on a vertical line, and not zero
    sum->zero = false;
    difnum(&(a->x),&(b->x),&xamxb);
    difnum(&(a->y),&(b->y),&yamyb);
    inversenum(&xamxb,&invxamxb);
    timesnum(&yamyb,&invxamxb,&s);
    timesnum(&s,&s,&s2);
    difnum(&s2,&(a->x),&s2mxa);
    difnum(&s2mxa,&(b->x),&(sum->x));
    difnum(&(a->x),&(sum->x),&xamxs);
    timesnum(&xamxs,&s,&stxamxs);
    difnum(&stxamxs,&(a->y),&(sum->y));  
    }
  /* we put the original modulo parameter back */
  for (i = 0 ; i < nofdigits ; i++)
    {
    modulo.digit[i] = oldmodulo.digit[i];
    }  
  } // end ellipsum
 

Finally, we use the square-and-multiply variant of the modular power, to implement the multiplication of an elliptic point with a number.  The number is transformed in a binary representation, and this binary representation is used to double ("square") the actual point, or to double the actual point, and add the original point to it ("multiply").

void elliptictimes(ellipcurve* mycurve, elliptic* x, number* y, elliptic* xy)
/* number x to the power y */
  {
  int i;
  int j;
  elliptic p;
  elliptic s;
  number zero;
  number half;
  number full;
  int uneven;
  int bitsequence[basisbits * nofdigits];
  int bitcounter;
  xy->zero = true;
  for (i = 0 ; i < nofdigits ; i++)
    {
    zero.digit[i] = 0;
    full.digit[i] = y->digit[i];
    }
  bitcounter = 0;
  /* we find the binary representation of the power */
  do
    {
    divide2(&full,&half,&uneven);
    bitsequence[bitcounter] = uneven;
    for (i = 0 ; i < nofdigits ; i++)
      {
      full.digit[i] = half.digit[i];
      }
    bitcounter++;
    }
  while (isbigger(&half,&zero) > 0);
  /* the square-and-multiply algorithm */
  for (i = bitcounter - 1 ; i >= 0 ; i--)
    {
    ellipsum(mycurve,xy,xy,&s);  // we "square"
    ellipsum(mycurve,&s,x,&p);   // we "multiply"
    if (bitsequence[i] == 1)
      {
      ellipcopy(&p,xy);
      }
    else
      {
      ellipcopy(&s,xy);
      }
    }
  }  // end elliptictimes
 

Before even being able to add two points or together, or to multiply a point with a number, we have to have a point.  Hasse's theorem indicates that the number of points in an elliptic group (the order N of the elliptic group) over a modular field with prime number p is not very far from p (the deviation of N from p is of the order of the square root of p).  However, we are working in the "plane" of p x p points, which means that a random x and a random y won't belong to the curve.

In order to find points on the curve, we can try an x-value.  The equation of the curve then gives us y2 ; however, that doesn't give us y - it doesn't even mean that a y exists such that the value squared gives us "y2".  If no such value exists, it means that the chosen x-value is not the x-value of any point on the curve.  Happily, it is not too difficult to test whether a given value for "y2" is, in fact, the square of an existing number in the modular field: the Legendre symbol.  If we want to know whether a given number x in a modular prime field is the square of another number (one calls x then a "quadratic residue"), one calculates x(p-1)/2.  If the result of this equals 1, then x is a quadratic residue ; if it equals -1, then x is not a quadratic residue.  The special case of 0 gives us 0 of course.  This is easily understood by considering Fermat's little theorem: if x = a2 then x(p-1)/2 = a(p-1) = 1.

The piece of code that does this calculation is given next:

int legendresymbol(number* prime,number* a)
  /* implementation of the Legendre symbol that verifies whether a
     is a square residue in the modular field with prime */
  {
  unsigned int i;
  number oldmodulo;
  number zero;
  number one;
  number mone;
  number pm1d2;
  number blegendre;
  int rr;
  int result;

  // we apply the prime number
  for (i = 0 ; i < nofdigits ; i++)
    {
    oldmodulo.digit[i] = modulo.digit[i];
    modulo.digit[i] = prime->digit[i];
    }
  makebig(0,&zero);
  makebig(1,&one);
  difnum(&zero,&one,&mone);
  divide2(&mone,&pm1d2,&rr);
  result = 15;
  if (rr != 0)
    {
    printf("legendresymbol: failure !\n");
    result = 2;
    }
    if (result == 15)
    {
    power(a,&pm1d2,&blegendre);
    if (isbigger(&blegendre,&one) == 0)
      {
      result = 1;
      }
    else if(isbigger(&blegendre,&zero) == 0)
      {
      result = 0;
      }
    else if (isbigger(&blegendre,&mone) == 0)
      {
      result = -1;
      }
    else
      {
      printf("legendresymbol: internal inconsistency!\n");
      result = 2;
      }
    }
  /* we put the original modulo parameter back */
  for (i = 0 ; i < nofdigits ; i++)
    {
    modulo.digit[i] = oldmodulo.digit[i];
    }  
  return result;
  }

Note that this function indicates whether a number is the square of another number, but that doesn't indicate us what that number ("the square root") actually is.  Finding the actual square roots if they exist, is somewhat more sophisticated, but there is a known algorithm, Cipolla's algorithm, that allows one to find them.

The algorithm needs a number a, such that a2 - x is not a quadratic residue.  Although there is no technique to construct such a number a, they are pretty abundant, and hence, a simple trial-and-error loop allows us to find such an a.  Once we have such an a, we consider the non-number sqrt(a2 - x) = w and we obtain a kind of "algebraic extension" when we include this non-number, in the same way as i = sqrt(-1) as a non-real-number gave rise to the complex numbers: the couple (u,v) over the modular prime field (like the couples of real numbers gave the complex numbers) is considered u + w.v ; but we keep the fact that w2 = a2 - x.

Capolla's algorithm simply consists in calculating (a + w)(p+1)/2.  This number doesn't contain any "w" component, and is claimed to be one of the square roots of x.  The other square root is of course minus this number.  In order to calculate this power in this extended algebra, we use a "square-and-multiply" technique, but we have to implement the multiplications and squares in the extended algebra:

  • (u,v)2 = (u2 + v2 w2 ; 2uv)
  • (u,v)(a+w) = (v w2 + a u ; u + a v)

Here is the implementation:

void cipolla(number* prime, number* square, number* root1, number* root2)
  {
  int i;
  unsigned int j;
  number oldmodulo;
  number a;
  number a2;
  number a2msquare;
  number zero;
  number one;
  number two;
  number pm1;
  number pm1d2;
  int rr;
  number pp1d2;
  number half;
  number full;
  int uneven;
  int bitsequence[basisbits * nofdigits];
  int bitcounter;
  number root1u;
  number root1v;
  number root1u2;
  number root1v2;
  number root1v2w2;
  number root1uv;
  number root1svw2;
  number root1sua;
  number root1sva;
  number su;
  number sv;
  number pu;
  number pv;
 
  // first, we check whether "square" is a quadratic residue
  if (legendresymbol(prime,square) != 1)
    {
    printf("cipolla: Trying to find a square root of a square non-residue !\n");
    return;
    }
  // we apply the prime number
  for (i = 0 ; i < nofdigits ; i++)
    {
    oldmodulo.digit[i] = modulo.digit[i];
    modulo.digit[i] = prime->digit[i];
    }  
  // finding a non-square residue a^2 - n
  do
    {
    randomgenmod(prime,&a);
    timesnum(&a,&a,&a2);
    difnum(&a2,square,&a2msquare);
    }
  while (legendresymbol(prime,&a2msquare) != -1);
  makebig(0,&zero);
  makebig(1,&one);
  makebig(2,&two);
  difnum(&zero,&one,&pm1);
  divide2(&pm1,&pm1d2,&rr);
  sumnum(&pm1d2,&one,&pp1d2);

  // preparing the square and multiply algorithm in F_p x F_p
  for (i = 0 ; i < nofdigits ; i++)
    {
    full.digit[i] = pp1d2.digit[i];
    }
  bitcounter = 0;
  /* we find the binary representation of the power */
  do
    {
    divide2(&full,&half,&uneven);
    bitsequence[bitcounter] = uneven;
    for (i = 0 ; i < nofdigits ; i++)
      {
      full.digit[i] = half.digit[i];
      }
    bitcounter++;
    }
  while (isbigger(&half,&zero) > 0);
  makebig(1,&root1u);
  makebig(0,&root1v);
 
  /* we will now find (a + 1 sqrt(a^2 - square))^(pp1d2)
  /* the square-and-multiply algorithm */
  for (i = bitcounter - 1 ; i >= 0 ; i--)
    {
    // we calculate root1^2 which is (u^2 + v^2*w^2, 2*u*v)
    timesnum(&root1u,&root1u,&root1u2);
    timesnum(&root1v,&root1v,&root1v2);
    timesnum(&root1v2,&a2msquare,&root1v2w2);
    sumnum(&root1u2,&root1v2w2,&su);
    timesnum(&root1u,&root1v,&root1uv);
    timesnum(&two,&root1uv,&sv);
    // we calculate s*(a,1) which is (sv*w^2 + su*a,su + a*sv)
    timesnum(&sv,&a2msquare,&root1svw2);
    timesnum(&a,&su,&root1sua);
    sumnum(&root1svw2,&root1sua,&pu);
    timesnum(&a,&sv,&root1sva);
    sumnum(&root1sva,&su,&pv);
    // we pick the right one for this bit
    if (bitsequence[i] == 1)
      {
      for (j = 0 ; j < nofdigits ; j++)
        {
        root1u.digit[j] = pu.digit[j];
        root1v.digit[j] = pv.digit[j];
        }
      }
    else
      {
      for (j = 0 ; j < nofdigits ; j++)
        {
        root1u.digit[j] = su.digit[j];
        root1v.digit[j] = sv.digit[j];
        }
      }
    }
  /* now we test whether the v-part of the root has become 0 */
  if (isbigger(&root1v,&zero) != 0)
    {
    printf("cipolla: inconsistency: v-part non-zero !\n");
    }

  /* root1u is one of the roots, -root1u is the other one */
  sumnum(&zero,&root1u,root1);
  difnum(&zero,&root1u,root2);
  /* we put the original modulo parameter back */
  for (i = 0 ; i < nofdigits ; i++)
    {
    modulo.digit[i] = oldmodulo.digit[i];
    }
  }  // end cipolla

One should only call this function if the square roots exist, by a test on the Legendre symbol.

Now that we have a way to take square roots in the modular prime field, we can finally find the possible y that go with a given x on an elliptic curve:

int ellipyfromx(ellipcurve* mycurve,number* x,number* y1,number* y2)
  /* calculates the two possible y-values for a given x-value on the curve
     if ever that x-value is possible.  The integer return value gives us the number
     of y-values (0,1 or 2) that are returned */
  {
  int i;
  number oldmodulo;
  number x2;
  number x3;
  number ax;
  number axpb;
  number yp2;
  int leg;
  int result;
 
  // we apply the prime number
  for (i = 0 ; i < nofdigits ; i++)
    {
    oldmodulo.digit[i] = modulo.digit[i];
    modulo.digit[i] = mycurve->prime.digit[i];
    }  
  // we calculate the y^2 corresponding to the x for the given curve
  timesnum(x,x,&x2);
  timesnum(&x2,x,&x3);
  timesnum(&(mycurve->a),x,&ax);
  sumnum(&ax,&(mycurve->b),&axpb);
  sumnum(&x3,&axpb,&yp2);
  leg = legendresymbol(&(mycurve->prime),&yp2);
  if (leg == -1)
    {
    result = 0;
    }
  else if (leg == 0)
    {
    result = 1;
    makebig(0,y1);
    makebig(0,y2);
    }
  else if (leg == 1)
    {
    cipolla(&(mycurve->prime),&yp2,y1,y2);
    result = 2;
    }
  else
    {
    result = -1;
    printf("ellipyfromx: inconsistency ! \n");
    }
    

  /* we put the original modulo parameter back */
  for (i = 0 ; i < nofdigits ; i++)
    {
    modulo.digit[i] = oldmodulo.digit[i];
    }
  return result;
  }  // end ellipyfromx

The code speaks for itself: we calculate u = x3 + ax + b, we verify if u is a square residue, and if it is, the two square roots of u are the corresponding possible y-values of the elliptic point with given x ; if not, x is not a good x-value of an elliptic point.  There is the odd case where u = 0, in which case, we only return y = 0.

We can finally present a toy elliptic curve generator.  We have to put a big warning here: this generator doesn't check on any of the weak curves, or on small cyclic groups.  This toy generator only serves to generate consistent, but maybe cryptographically weak, elliptic groups.  We only use it to illustrate cryptographic algorithms, and to help do automatic random consistency checks on the code.

void toycurve(ellipcurve* mycurve,int nbits)
  /* inefficient algorithm to generate small elliptic curves */
  /* impossible to use for real crypto */
  {

  number oldmodulo;
  number zero;
  number one;
  number two;
  number nbitsbig;
  number allbits;
  number ranlim;
  int i;
  int yok;
  number dummy;

  /* we put the current modulus away and set the modulus to its maximum */
  for (i = 0 ; i < nofdigits ; i++)
    {
    oldmodulo.digit[i] = modulo.digit[i];
    modulo.digit[i] = basis - 1;
    zero.digit[i] = 0;
    one.digit[i] = 0;
    two.digit[i] = 0;
    }
  one.digit[0] = 1;
  two.digit[0] = 2;
  /* we calculate the highest possible number in agreement with our number of bits (minus 1) */
  makebig(nbits-1,&nbitsbig);
  power(&two,&nbitsbig,&allbits);
  /* we just find a prime number of the order of n bits */
  do
    {
    randomgenmod(&allbits,&ranlim);
    findnextprimerabin(&ranlim,&(mycurve->prime),10);
    }
  while (isbigger(&(mycurve->prime),&allbits) == 1);

  /* we generate a random curve (random a and b) */
  randomgenmod(&(mycurve->prime),&(mycurve->a));
  randomgenmod(&(mycurve->prime),&(mycurve->b));
 
  /* we now have to find a random point on the curve.
     For this, we have to generate a random x value, and
     verify whether the corresponding y^2 = x^3 + ax + b
     is a square residue */
  mycurve->generator.zero = false;
  do
    {
    randomgenmod(&(mycurve->prime),&(mycurve->generator.x));
    yok = ellipyfromx(mycurve,&(mycurve->generator.x),&(mycurve->generator.y),&dummy);
    }
  while (yok < 2);

  /* we put the original modulo parameter back */
  for (i = 0 ; i < nofdigits ; i++)
    {
    modulo.digit[i] = oldmodulo.digit[i];
    }
  }  // end toycurve

Again, toycurve generates a prime field, generates a random elliptic curve on that prime field, and then picks a random generator on that curve, which will establish a cyclic group. No single check is done to eliminate all the different weak cases, or the small cyclic groups, or the cyclic groups with small subgroups, which will all lead to a breakable elliptic curve cryptography.