/*
  Author:  Pate Williams (c) 1997

  Example
  "6.6 As before, suppose p = 467, alpha = 4,
  a = 101 and beta = 449. Suppose the message
  x = 286 is signed with the (bogus) signature
  y = 83, and Bob wants to convince Alice that
  the signature is invalid."
  -Douglas R. Stinson-
  See "Cryptogrphy: Theory and Practice" by
  Douglas R. Stinson pages 222-223.
*/

#include <stdio.h>

long exp_mod(long x, long b, long n)
/* returns x ^ b mod n */
{
  long a = 1l, s = x;

  while (b != 0) {
    if (b & 1l) a = (a * s) % n;
    b >>= 1;
    if (b != 0) s = (s * s) % n;
  }
  if (a < 0) a += n;
  return a;
}

long Extended_Euclidean(long b, long n)
{
  long b0 = b, n0 = n, t = 1, t0 = 0, temp, q, r;

  q = n0 / b0;
  r = n0 - q * b0;
  while (r > 0) {
    temp = t0 - q * t;
    if (temp >= 0) temp = temp % n;
    else temp = n - (- temp % n);
    t0 = t;
    t = temp;
    n0 = b0;
    b0 = r;
    q = n0 / b0;
    r = n0 - q * b0;
  }
  if (b0 != 1) return 0;
  else return t % n;
}

int main(void)
{
  long a = 101, alpha = 4, beta = 449, e1 = 45;
  long e2 = 237, f1 = 125, f2 = 9, i, j, p = 467;
  long q, x = 286, y = 83, c, d, C, D, r, s, t;

  q = (p - 1) >> 1;
  printf("a = %ld\n", a);
  printf("alpha = %ld\n", alpha);
  printf("beta = %ld\n", beta);
  printf("e1 = %ld\n", e1);
  printf("e2 = %ld\n", e2);
  printf("f1 = %ld\n", f1);
  printf("f2 = %ld\n", f2);
  printf("p = %ld\n", p);
  printf("q = %ld\n", q);
  printf("x = %ld\n", x);
  printf("y = %ld\n", y);
  i = Extended_Euclidean(a, q);
  c = (exp_mod(y, e1, p) * exp_mod(beta, e2, p)) % p;
  d = exp_mod(c, i, p);
  printf("Alice's challenge c = %ld\n", c);
  printf("Bob's response d = %ld\n", d);
  if (d != (exp_mod(x, e1, p) * exp_mod(alpha, e2, p)) % p)
    printf("d != x ^ e1 * alpha ^ e2 mod p\n");
  else
    printf("d == x ^ e1 * alpha ^ e2 mod p\n");
  C = (exp_mod(y, f1, p) * exp_mod(beta, f2, p)) % p;
  D = exp_mod(C, i, p);
  printf("Alice's challenge C = %ld\n", C);
  printf("Bob's response D = %ld\n", D);
  if (D != (exp_mod(x, f1, p) * exp_mod(alpha, f2, p)) % p)
    printf("D != x ^ f1 * alpha ^ f2 mod p\n");
  else
    printf("D == x ^ f1 * alpha ^ f2 mod p\n");
  i = q - e2;
  if (i < 0) i += q;
  j = q - f2;
  if (j < 0) j += q;
  r = (d * exp_mod(alpha, i, p)) % p;
  s = exp_mod(r, f1, p);
  r = (D * exp_mod(alpha, j, p)) % p;
  t = exp_mod(r, e1, p);
  if (s == t)
    printf("Alice concludes y is a forgery\n");
  else
    printf("Alice does not conclude y is a forgery\n");
  return 0;
}