/*

   Scott's AKE Client/Server testbed

   See http://eprint.iacr.org/2002/164

   Compile as 
   cl /O2 /GX /DZZNS=5 ake6mntt.cpp zzn6.cpp ecn3.cpp zzn3.cpp big.cpp zzn.cpp ecn.cpp miracl.lib
   using COMBA build

   MNT Curve - Tate pairing

   The required file mnt.ecs is created from a curve generated by the mnt 
   utility, and created by the cm utility. For convenience the value of 
   (p^2-p+1)/q and the 6th root of unity (cnr^(p-1)/6) have been manually 
   calculated and appended to this file (replacing the x,y values in the
   original .ecs file)

   NOTE: Irreducible polynomial MUST be of the form x^6+CNR. This excludes many of the curves
   found using the mnt utility!

   Use the irred utility

   Modified to prevent sub-group confinement attack

   NOTE: Key exchange bandwidth could be reduced further using ideas from
         "Doing more with Fewer Bits", Brouwer, Pellikaan & Verheul, Asiacrypt 
         '99

   Speeded up using ideas from 
   "Efficient Computation of Tate Pairing in Projective Coordinate over General Characteristic Fields" 
   by Sanjit Chatterjee1, Palash Sarkar1 and Rana Barua1 

*/

#include <iostream>
#include <fstream>
#include <string.h>
#include "ecn.h"
#include <ctime>
#include "ecn3.h"
#include "zzn6.h"

// fix a couple of things for this particular curve
// cofactor - number of points on curve=CF.q
// Cubic non-residue mod p

#define CF 2  
#define CNR 2 // irreducible is x^6-2

using namespace std;

Miracl precision(5,0); 

#ifdef MR_COUNT_OPS
extern "C"
{
    int fpc,fpa,fpx,fpm2,fpi2;
}
#endif

// Using SHA-1 as basic hash algorithm

#define HASH_LEN 20

//
// Define one or the other of these
//
// Which is faster depends on the I/M ratio - See imratio.c
// Roughly if I/M ratio > 16 use PROJECTIVE, otherwise use AFFINE
//

#ifdef MR_AFFINE_ONLY
    #define AFFINE
#else
    #define PROJECTIVE
#endif

//
// Tate Pairing Code
//
// Extract ECn point in internal ZZn format
//

void extract(ECn& A,ZZn& x,ZZn& y)
{ 
    x=(A.get_point())->X;
    y=(A.get_point())->Y;
}

#ifdef PROJECTIVE
void extract(ECn& A,ZZn& x,ZZn& y,ZZn& z)
{ 
    big t;
    x=(A.get_point())->X;
    y=(A.get_point())->Y;
    t=(A.get_point())->Z;
    if (A.get_status()!=MR_EPOINT_GENERAL) z=1;
    else                                   z=t;
}
#endif

//
// Line from A to destination C. Let A=(x,y)
// Line Y-slope.X-c=0, through A, so intercept c=y-slope.x
// Line Y-slope.X-y+slope.x = (Y-y)-slope.(X-x) = 0
// Now evaluate at Q -> return (Qy-y)-slope.(Qx-x)
//

ZZn6 line(ECn& A,ECn& C,ECn& B,int type,ZZn& slope,ZZn& ex1,ZZn& ex2,ZZn3& Qx,ZZn3& Qy)
{
     ZZn6 w;
#ifdef AFFINE
     ZZn3 nn=Qx;
     ZZn x,y;
     extract(A,x,y);
     nn-=x; 
     nn*=slope;   
     nn+=y;
     w.set(-nn,Qy);
#endif
#ifdef PROJECTIVE
     if (type==MR_ADD)
     {
        ZZn x2,y2,x3,z3;
        extract(B,x2,y2);
        extract(C,x3,x3,z3);
        w.set(slope*(Qx-x2)+z3*y2,-z3*Qy);
     } 
     if (type==MR_DOUBLE)
     {
        ZZn x,y,x3,z3;
        extract(A,x,y);
        extract(C,x3,x3,z3);
        w.set((slope*ex2)*Qx-slope*x+ex1,-(z3*ex2)*Qy);
     }
/*     extract(A,x,y,z);
     x*=z; t=z; z*=z; z*=t;
     x*=slope; t=slope*z;
     nn*=t; nn-=x; t=z;
     extract(C,x,x,z);
     nn+=(z*y); t*=z;
     w.set(nn,-Qy*t);
*/
#endif
     return w;
}

//
// Add A=A+B  (or A=A+A) 
// Return line function value
//

ZZn6 g(ECn& A,ECn& B,ZZn3& Qx,ZZn3& Qy)
{
    ZZn  lam,extra1,extra2;
    int type;
    ZZn6 u;
    big ptr,ex1,ex2;
    ECn P=A;

// Evaluate line from A
    type=A.add(B,&ptr,&ex1,&ex2);
    if (!type) return (ZZn6)1;
    lam=ptr;
    extra1=ex1;
    extra2=ex2;
    return line(P,A,B,type,lam,extra1,extra2,Qx,Qy);
}

//
// Tate Pairing - note denominator elimination has been applied
//
// P is a point of order q. Q(x,y) is a point of order m.q. 
// Note that P is a point on the curve over Fp, Q(x,y) a point on the 
// twisted curve over the extension field Fp^3
//

BOOL fast_tate_pairing(ECn& P,ZZn3& Qx,ZZn3& Qy,Big& q,Big &cf,ZZn6& res)
{ 
    int i,j,n,nb,nbw,nzs;
    ECn A,P2,t[8];
    ZZn6 w,hc,z2n,zn[8];

    res=zn[0]=1;  

    t[0]=P2=A=P;
	z2n=g(P2,P2,Qx,Qy);

    normalise(P2);
//
// Build windowing table
//

    for (i=1;i<8;i++)
    {           
        hc=g(A,P2,Qx,Qy);      
        t[i]=A;         
        zn[i]=z2n*zn[i-1]*hc;   
    }
  
    multi_norm(8,t);  // make t points Affine

    A=P;    // reset A
    nb=bits(q);

    for (i=nb-2;i>=0;i-=(nbw+nzs))
    {
        n=window(q,i,&nbw,&nzs,4);  // standard MIRACL windowing

        for (j=0;j<nbw;j++)
        {
            res*=res;   
            res*=g(A,A,Qx,Qy);
        }

        if (n>0)
        {
            res*=zn[n/2]; 
            res*=g(A,t[n/2],Qx,Qy);
        }  

        for (j=0;j<nzs;j++) 
        {
            res*=res; 
            res*=g(A,A,Qx,Qy); 
        } 
	
        if (res.iszero()) return FALSE;  
    }
#ifdef MR_COUNT_OPS
printf("After Miller fpc= %d fpa= %d fpx= %d\n",fpc,fpa,fpx);
#endif
  //  if (!A.iszero() || res.iszero()) return FALSE;

    w=res;   
    w.powq();
    res*=w;                        // ^(p+1)

    w=res;
    w.powq(); w.powq(); w.powq();
    res=w/res;                     // ^(p^3-1)

// res=pow(res,cf);               // ^(p*p-p+1)/q

// exploit the clever "trick" for a half-length exponentiation!
// cout << "Final Exponentiation" << endl;
    res.mark_as_unitary();

    w=res;
    res.powq(); res*=res; // res*=res;  // res=pow(res,CF);
 
    if (cf<0) res/=powu(w,-cf);
    else res*=powu(w,cf);
//cout << "res= " << res << endl;
//exit(0);


/*
cout << "res*res= " << res*res << endl;

res.make_unitary();

cout << "res*res= " << res*res << endl;

ZZn3 r0,r1;

r0=real(res);
r1=imaginary(res);

cout << "real= " << 2*tx(r1*r1)+(ZZn3)1 << endl;
cout << "Imaginary= " << (r0+r1)*(r0+r1)-tx(r1*r1)-(ZZn3)1-r1*r1 << endl;


exit(0);


w=res;

r=2*real(res);

ZZn3 rp,ra;

res.powq();
rp=2*real(res);

res/=w;
ra=2*real(res);

Big a=rand(q);

cout << "powl(r,a)= " << powl(r,a) << endl;

Big a0,a1;
a0=a%T; a1=a/T;

cout << "pow(r,a)= " << powl(rp,a1,r,a0,ra) << endl;

    exit(0);
*/

    if (res==(ZZn6)1) return FALSE;
    return TRUE;            
}

//
// ecap(.) function
//

BOOL ecap(ECn& P,ECn3& Q,Big& order,Big& cf,ZZn6& res)
{
    BOOL Ok;
    ECn PP=P;
    ZZn3 Qx,Qy;
    int qnr=get_mip()->cnr;

    normalise(PP);
    Q.get(Qx,Qy);

// untwist    
    Qx=Qx/qnr;
    Qy=tx(Qy);
    Qy=Qy/(qnr*qnr);

#ifdef MR_COUNT_OPS
fpc=fpa=fpx=0;
#endif  

    Ok=fast_tate_pairing(PP,Qx,Qy,order,cf,res);

#ifdef MR_COUNT_OPS
printf("After tate fpc= %d fpa= %d fpx= %d\n",fpc,fpa,fpx);
fpa=fpc=fpx=0;
#endif

    if (Ok) return TRUE;
    return FALSE;
}

//
// Hash functions
// 

Big H1(char *string)
{ // Hash a zero-terminated string to a number < modulus
    Big h,p;
    char s[HASH_LEN];
    int i,j; 
    sha sh;

    shs_init(&sh);

    for (i=0;;i++)
    {
        if (string[i]==0) break;
        shs_process(&sh,string[i]);
    }
    shs_hash(&sh,s);
    p=get_modulus();
    h=1; j=0; i=1;
    forever
    {
        h*=256; 
        if (j==HASH_LEN)  {h+=i++; j=0;}
        else         h+=s[j++];
        if (h>=p) break;
    }
    h%=p;
    return h;
}

Big H2(ZZn3 x)
{ // Hash an Fp3 to a big number
    sha sh;
    ZZn u,v,w;
    Big a,h,p,xx[3];
    char s[HASH_LEN];
    int i,j,m;

    shs_init(&sh);
    x.get(u,v,w);
    xx[0]=u; xx[1]=v; xx[2]=w;
    for (i=0;i<3;i++)
    {
        a=xx[i];
        while (a>0)
        {
            m=a%256;
            shs_process(&sh,m);
            a/=256;
        }
    }
    shs_hash(&sh,s);
    h=from_binary(HASH_LEN,s);
    return h;
}

// Hash and map a Server Identity to a curve point E_(Fp3)

ECn3 hash_and_map3(char *ID)
{
    int i;
    ECn3 S;
    ZZn3 X;
 
    Big x0=H1(ID);
    forever
    {
        x0+=1;
        X.set2((ZZn)x0);
        if (!S.set(X)) continue;

        break;
    }

//    cout << "S= " << S << endl;
    return S;
}     

// Hash and map a Client Identity to a curve point E_(Fp) of order q

ECn hash_and_map(char *ID)
{
    ECn Q;
    Big x0=H1(ID);

    while (!Q.set(x0,x0)) x0+=1;
    Q*=CF;
    return Q;
}

// Use Galbraith & Scott Homomorphism idea

ZZn3 mypow(ZZn6& res,Big &e,Big &T)
{
	ZZn6 w=res;
	ZZn3 ra,rp,r=real(res);
	Big e0,e1;
	e0=e%T; e1=e/T;
	w.powq(); rp=real(w); w/=res; ra=real(w);
// Use GLV method, and double exponentiation a la Lucas (see ZZn3.cpp)
	return powl(rp,e1,r,e0,ra);
}

int main()
{
    ifstream common("mnt.ecs");      // MNT elliptic curve parameters
    miracl* mip=&precision;
    ECn Alice,Bob,sA,sB;
    ECn3 B6,Server,sS;
    ZZn3 sp,ap,bp;
	ZZn6 res;
    Big a,b,s,ss,p,q,x,y,B,cf,cfp,t,sru,T;
    int i,bits,A;
    time_t seed;

    common >> bits;
    mip->IOBASE=16;
    common >> p;
    common >> A;
    common >> B >> q >> cf >> sru;

	T=p-q*CF;

    time(&seed);
    irand((long)seed);

#ifdef AFFINE
    ecurve(A,B,p,MR_AFFINE);
#endif
#ifdef PROJECTIVE
    ecurve(A,B,p,MR_PROJECTIVE);
#endif

    set_zzn3(CNR,sru);
    cfp=cf-CF*p;  // ~ (t-1)

    mip->IOBASE=16;
    mip->TWIST=MR_QUADRATIC;   // map Server to point on twisted curve E(Fp3)

    ss=rand(q);    // TA's super-secret 

    cout << "Mapping Server ID to point" << endl;
    Server=hash_and_map3((char *)"Server");
    sS=ss*Server;

    cout << "Mapping Alice & Bob ID's to points" << endl;
    Alice=hash_and_map((char *)"Alice");
    Bob=  hash_and_map((char *)"Robert");

    cout << "Alice, Bob and the Server visit Trusted Authority" << endl; 


    sA=ss*Alice; 
    sB=ss*Bob; 

    cout << "Alice and Server Key Exchange" << endl;

    a=rand(q);   // Alice's random number
    s=rand(q);   // Server's random number

    if (!ecap(sA,Server,q,cfp,res)) cout << "Trouble" << endl;

    if (powl(real(res),q)!=(ZZn3)1)
    {
        cout << "Wrong group order - aborting" << endl;
        exit(0);
    }
//    ap=powl(real(res),a);
	ap=mypow(res,a,T);

//for (i=0;i<10000;i++)
    if (!ecap(Alice,sS,q,cfp,res)) cout << "Trouble" << endl;
    if (powl(real(res),q)!=(ZZn3)1)
    {
        cout << "Wrong group order - aborting" << endl;
        exit(0);
    }
//    sp=powl(real(res),s);
	sp=mypow(res,s,T);

    cout << "Alice  Key= " << H2(powl(sp,a)) << endl;
    cout << "Server Key= " << H2(powl(ap,s)) << endl;

    cout << "Bob and Server Key Exchange" << endl;

    b=rand(q);   // Bob's random number
    s=rand(q);   // Server's random number

    if (!ecap(sB,Server,q,cfp,res)) cout << "Trouble" << endl;
    if (powl(real(res),q)!=(ZZn3)1)
    {
        cout << "Wrong group order - aborting" << endl;
        exit(0);
    }
    bp=powl(real(res),b);

    if (!ecap(Bob,sS,q,cfp,res)) cout << "Trouble" << endl;
    if (powl(real(res),q)!=(ZZn3)1)
    {
        cout << "Wrong group order - aborting" << endl;
        exit(0);
    }
    sp=powl(real(res),s);

    cout << "Bob's  Key= " << H2(powl(sp,b)) << endl;
    cout << "Server Key= " << H2(powl(bp,s)) << endl;

    return 0;
}