mrcurve.c

大数运算库

/*
 *   MIRACL elliptic curve routines
 *   mrcurve.c
 *
 *   Assumes Weierstrass equation y^2 = x^3 + Ax + B
 *   See IEEE P1363 Draft Standard 
 *
 *   Uses Montgomery's representation internally
 *
 *   Works particularly well with fixed length Comba multiplier
 *   e.g. #define MR_COMBA 5 for 5x32 = 160 bit modulus
 *        on 32-bit computer
 *
 *   Copyright (c) 1997-2001 Shamus Software Ltd.
 */

#include <stdio.h>
#include <miracl.h>

#define mr_abs(x)  ((x)<0? (-(x)) : (x))

/* initialise elliptic curve */

void ecurve_init(_MIPD_ big a,big b,big p,int type)
{ /* Initialize the active ecurve    *
   * Asize indicate size of A        *
   * Bsize indicate size of B        */
    int as;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return;

    MR_IN(93)

    mr_mip->SS=FALSE;       /* no special support for super-singular curves */ 
    prepare_monty(_MIPP_ p);

    mr_mip->Asize=size(a);
    if (mr_abs(mr_mip->Asize)==MR_TOOBIG)
    {
        if (mr_mip->Asize>=0)
        { /* big positive number - check it isn't minus something small */
           copy(a,mr_mip->w1);
           divide(_MIPP_ mr_mip->w1,p,p);
           subtract(_MIPP_ p,mr_mip->w1,mr_mip->w1);
           as=size(mr_mip->w1);
           if (as<MR_TOOBIG) mr_mip->Asize=-as;
           else
           {
               if (mr_mip->A==NULL) mr_mip->A=mirvar(_MIPP_ 0);
               nres(_MIPP_ a,mr_mip->A);
           }
        }
        else
        {
               if (mr_mip->A==NULL) mr_mip->A=mirvar(_MIPP_ 0);
               nres(_MIPP_ a,mr_mip->A);
        }
    }
    mr_mip->Bsize=size(b);
    if (mr_abs(mr_mip->Bsize)==MR_TOOBIG)
    {
        if (mr_mip->B==NULL) mr_mip->B=mirvar(_MIPP_ 0);
        nres(_MIPP_ b,mr_mip->B);

    }
    mr_mip->coord=type;
    MR_OUT
    return;
}

epoint* epoint_init(_MIPDO_ )
{ /* initialise epoint to point at infinity. */
    epoint *p;
    char *ptr;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return NULL;

    MR_IN(96)

/* Create space for whole structure in one heap access */ 

    if (mr_mip->coord!=MR_AFFINE)
        p=(epoint *)mr_alloc(_MIPP_ sizeof(epoint)+3*mr_mip->size,1);
    else
        p=(epoint *)mr_alloc(_MIPP_ sizeof(epoint)+2*mr_mip->size,1);
        
    ptr=(char *)p+sizeof(epoint);
    p->X=mirvar_mem(_MIPP_ ptr,0);
    p->Y=mirvar_mem(_MIPP_ ptr,1);
    if (mr_mip->coord!=MR_AFFINE) p->Z=mirvar_mem(_MIPP_ ptr,2);
    else p->Z=NULL;

    p->marker=MR_EPOINT_INFINITY;

    MR_OUT

    return p;
}

void epoint_free(epoint *p)
{ /* clean up point */
    zero(p->X);
    zero(p->Y);
    if (p->marker==MR_EPOINT_GENERAL) zero(p->Z);
    mr_free(p);
}        

BOOL epoint_set(_MIPD_ big x,big y,int cb,epoint *p)
{ /* initialise a point on active ecurve            *
   * if x or y == NULL, set to point at infinity    *
   * if x==y, a y co-ordinate is calculated - if    *
   * possible - and cb suggests LSB 0/1  of y       *
   * (which "decompresses" y). Otherwise, check     *
   * validity of given (x,y) point, ignoring cb.    *
   * Returns TRUE for valid point, otherwise FALSE. */
  
    BOOL valid;

#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return FALSE;

    MR_IN(97)

    if (x==NULL || y==NULL)
    {
        convert(_MIPP_ 1,mr_mip->w1);
        nres(_MIPP_ mr_mip->w1,p->X);
        nres(_MIPP_ mr_mip->w1,p->Y);
        p->marker=MR_EPOINT_INFINITY;
        MR_OUT
        return TRUE;
    }

/* find x^3+Ax+B */
    nres(_MIPP_ x,p->X);

    nres_modmult(_MIPP_ p->X,p->X,mr_mip->w3);

    nres_modmult(_MIPP_ mr_mip->w3,p->X,mr_mip->w3);
    if (mr_abs(mr_mip->Asize)==MR_TOOBIG)
        nres_modmult(_MIPP_ p->X,mr_mip->A,mr_mip->w1);
    else
        nres_premult(_MIPP_ p->X,mr_mip->Asize,mr_mip->w1);
    nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w1,mr_mip->w3);
    if (mr_abs(mr_mip->Bsize)==MR_TOOBIG)
        nres_modadd(_MIPP_ mr_mip->w3,mr_mip->B,mr_mip->w3);
    else
    {
        convert(_MIPP_ mr_mip->Bsize,mr_mip->w1);
        nres(_MIPP_ mr_mip->w1,mr_mip->w1);
        nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w1,mr_mip->w3);
    }
    valid=FALSE;
    if (x!=y)
    { /* compare with y^2 */
        nres(_MIPP_ y,p->Y);
        nres_modmult(_MIPP_ p->Y,p->Y,mr_mip->w1);
        
        if (compare(mr_mip->w1,mr_mip->w3)==0) valid=TRUE;
    }
    else
    { /* no y supplied - calculate one. Find square root */

        valid=nres_sqroot(_MIPP_ mr_mip->w3,p->Y);

    /* check LSB - have we got the right root? */
        redc(_MIPP_ p->Y,mr_mip->w1);
        if (remain(_MIPP_ mr_mip->w1,2)!=cb) 
            mr_psub(_MIPP_ mr_mip->modulus,p->Y,p->Y);

    } 
    if (valid)
    {
        p->marker=MR_EPOINT_NORMALIZED;
        MR_OUT
        return TRUE;
    }

    MR_OUT
    return FALSE;
}

int epoint_get(_MIPD_ epoint* p,big x,big y)
{ /* Get point co-ordinates in affine, normal form       *
   * (converted from projective, Montgomery form)        *
   * if x==y, supplies x only. Return value is Least     *
   * Significant Bit of y (useful for point compression) */

    int lsb;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif

    if (p->marker==MR_EPOINT_INFINITY)
    {
        zero(x);
        zero(y);
        return 0;
    }
    if (mr_mip->ERNUM) return 0;

    MR_IN(98)

    if (!epoint_norm(_MIPP_ p)) 
    { /* not possible ! */
        MR_OUT
        return (-1);
    }

    redc(_MIPP_ p->X,x);
    redc(_MIPP_ p->Y,mr_mip->w1);

    if (x!=y) copy(mr_mip->w1,y);
    lsb=remain(_MIPP_ mr_mip->w1,2); 
    MR_OUT
    return lsb;
}

BOOL epoint_norm(_MIPD_ epoint *p)
{ /* normalise a point */
    
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->coord==MR_AFFINE) return TRUE;
    if (p->marker!=MR_EPOINT_GENERAL) return TRUE;

    if (mr_mip->ERNUM) return FALSE;

    MR_IN(117)

    convert(_MIPP_ 1,mr_mip->w8);
    nres(_MIPP_ mr_mip->w8,mr_mip->w8);  

    if (nres_moddiv(_MIPP_ mr_mip->w8,p->Z,mr_mip->w8)>1) /* 1/Z  */
    {
        MR_OUT
        return FALSE;
    }
    
    nres_modmult(_MIPP_ mr_mip->w8,mr_mip->w8,mr_mip->w1);/* 1/ZZ */
    nres_modmult(_MIPP_ p->X,mr_mip->w1,p->X);            /* X/ZZ */
    nres_modmult(_MIPP_ mr_mip->w1,mr_mip->w8,mr_mip->w1);/* 1/ZZZ */
    nres_modmult(_MIPP_ p->Y,mr_mip->w1,p->Y);            /* Y/ZZZ */

    convert(_MIPP_ 1,mr_mip->w8);
    nres(_MIPP_ mr_mip->w8,p->Z);
   
    p->marker=MR_EPOINT_NORMALIZED;
    MR_OUT
    return TRUE;
}

void ecurve_multi_add(_MIPD_ int m,epoint **x,epoint**w)
{ /* adds m points together simultaneously, w[i]+=x[i] */
    int i,*flag;
    big *A,*B,*C;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return;    

    MR_IN(122)
    if (mr_mip->coord==MR_AFFINE)
    { /* this can be done faster */
        A=(big *)mr_alloc(_MIPP_ m,sizeof(big));
        B=(big *)mr_alloc(_MIPP_ m,sizeof(big));
        C=(big *)mr_alloc(_MIPP_ m,sizeof(big));
        flag=(int *)mr_alloc(_MIPP_ m,sizeof(int));

        convert(_MIPP_ 1,mr_mip->w3);
        nres(_MIPP_ mr_mip->w3,mr_mip->w3);    /* unity */
      
        for (i=0;i<m;i++)
        {
            A[i]=mirvar(_MIPP_ 0);
            B[i]=mirvar(_MIPP_ 0);
            C[i]=mirvar(_MIPP_ 0);
            flag[i]=0;
            if (compare(x[i]->X,w[i]->X)==0 && compare(x[i]->Y,w[i]->Y)==0) 
            { /* doubling */
                if (x[i]->marker==MR_EPOINT_INFINITY || size(x[i]->Y)==0)
                {
                    flag[i]=1;       /* result is infinity */
                    copy(mr_mip->w3,B[i]);
                    continue;    
                }
                nres_modmult(_MIPP_ x[i]->X,x[i]->X,A[i]);
                nres_premult(_MIPP_ A[i],3,A[i]);  /* 3*x^2 */
                if (mr_abs(mr_mip->Asize) == MR_TOOBIG)
                    nres_modadd(_MIPP_ A[i],mr_mip->A,A[i]);
                else
                {
                    convert(_MIPP_ mr_mip->Asize,mr_mip->w2);
                    nres(_MIPP_ mr_mip->w2,mr_mip->w2);
                    nres_modadd(_MIPP_ A[i],mr_mip->w2,A[i]);
                }                                       /* 3*x^2+A */
                nres_premult(_MIPP_ x[i]->Y,2,B[i]);
            }
            else
            {
                if (x[i]->marker==MR_EPOINT_INFINITY)
                {
                    flag[i]=2;              /* w[i] unchanged */
                    copy(mr_mip->w3,B[i]);
                    continue;
                }
                if (w[i]->marker==MR_EPOINT_INFINITY)
                {
                    flag[i]=3;              /* w[i] = x[i] */
                    copy(mr_mip->w3,B[i]);
                    continue;
                }
                nres_modsub(_MIPP_ x[i]->X,w[i]->X,B[i]);
                if (size(B[i])==0)
                { /* point at infinity */
                    flag[i]=1;       /* result is infinity */
                    copy(mr_mip->w3,B[i]);
                    continue;    
                }
                nres_modsub(_MIPP_ x[i]->Y,w[i]->Y,A[i]);
            }   
        }
        nres_multi_inverse(_MIPP_ m,B,C);  /* only one inversion needed */
        for (i=0;i<m;i++)
        {
            if (flag[i]==1)
            { /* point at infinity */
                epoint_set(_MIPP_ NULL,NULL,0,w[i]);
                continue;
            }
            if (flag[i]==2)
            {
                continue;
            }
            if (flag[i]==3)
            {
                epoint_copy(x[i],w[i]);
                continue;
            }
            nres_modmult(_MIPP_ A[i],C[i],mr_mip->w8);

            nres_modmult(_MIPP_ mr_mip->w8,mr_mip->w8,mr_mip->w2); /* m^2 */
            nres_modsub(_MIPP_ mr_mip->w2,x[i]->X,mr_mip->w1);
            nres_modsub(_MIPP_ mr_mip->w1,w[i]->X,mr_mip->w1);
       
            nres_modsub(_MIPP_ w[i]->X,mr_mip->w1,mr_mip->w2);
            nres_modmult(_MIPP_ mr_mip->w2,mr_mip->w8,mr_mip->w2);
            nres_modsub(_MIPP_ mr_mip->w2,w[i]->Y,w[i]->Y);
            copy(mr_mip->w1,w[i]->X);
            w[i]->marker=MR_EPOINT_GENERAL;

            mr_free(C[i]);
            mr_free(B[i]);
            mr_free(A[i]);
        }
        mr_free(flag);
        mr_free(C); mr_free(B); mr_free(A);
    }
    else
    { /* no speed-up */
        for (i=0;i<m;i++) ecurve_add(_MIPP_ x[i],w[i]);
    }
    MR_OUT
}

static void ecurve_double(_MIPD_ epoint *p)
{ /* double epoint on active ecurve */

#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (p->marker==MR_EPOINT_INFINITY) 
    { /* 2 times infinity == infinity ! */
        return;
    }

    if (mr_mip->coord==MR_AFFINE)
    {
        if (size(p->Y)==0) 
        { /* set to point at infinity */
            epoint_set(_MIPP_ NULL,NULL,0,p);
            return;
        }
        nres_modmult(_MIPP_ p->X,p->X,mr_mip->w8);    /* w8=x^2   */
        nres_premult(_MIPP_ mr_mip->w8,3,mr_mip->w8); /* w8=3*x^2 */
        if (mr_abs(mr_mip->Asize) == MR_TOOBIG)
            nres_modadd(_MIPP_ mr_mip->w8,mr_mip->A,mr_mip->w8);
        else
        {
            convert(_MIPP_ mr_mip->Asize,mr_mip->w2);
            nres(_MIPP_ mr_mip->w2,mr_mip->w2);
            nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w2,mr_mip->w8);
        }                                     /* w8=3*x^2+A */
        nres_premult(_MIPP_ p->Y,2,mr_mip->w6);      /* w6=2y */
        nres_moddiv(_MIPP_ mr_mip->w8,mr_mip->w6,mr_mip->w8);  

/* w8 is slope m on exit */

        nres_modmult(_MIPP_ mr_mip->w8,mr_mip->w8,mr_mip->w2); /* w2=m^2 */
        nres_premult(_MIPP_ p->X,2,mr_mip->w1);
        nres_modsub(_MIPP_ mr_mip->w2,mr_mip->w1,mr_mip->w1); /* w1=m^2-2x */
        
        nres_modsub(_MIPP_ p->X,mr_mip->w1,mr_mip->w2);
        nres_modmult(_MIPP_ mr_mip->w2,mr_mip->w8,mr_mip->w2);
        nres_modsub(_MIPP_ mr_mip->w2,p->Y,p->Y);
        copy(mr_mip->w1,p->X);
        
        return;
    }

    if (size(p->Y)==0)
    { /* set to point at infinity */
        epoint_set(_MIPP_ NULL,NULL,0,p);
        return;
    }
    convert(_MIPP_ 1,mr_mip->w1);
    if (p->marker==MR_EPOINT_NORMALIZED) nres(_MIPP_ mr_mip->w1,mr_mip->w2);
    else nres_modmult(_MIPP_ p->Z,p->Z,mr_mip->w2);
    if (mr_abs(mr_mip->Asize) < MR_TOOBIG)
    {
        if (mr_mip->Asize==(-3))
        { /* a is -3. Goody. */
            nres_modsub(_MIPP_ p->X,mr_mip->w2,mr_mip->w3);
            nres_modadd(_MIPP_ p->X,mr_mip->w2,mr_mip->w8);
            nres_modmult(_MIPP_ mr_mip->w3,mr_mip->w8,mr_mip->w3);
            nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w3,mr_mip->w8);
            nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w3,mr_mip->w8);
        }
        else
        { /* a is small */
            if (mr_mip->Asize!=0)
            { /* a is non zero! */
                nres_modmult(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w3);
                nres_premult(_MIPP_ mr_mip->w3,mr_mip->Asize,mr_mip->w3);
            }
            nres_modmult(_MIPP_ p->X,p->X,mr_mip->w1);
            nres_modadd(_MIPP_ mr_mip->w1,mr_mip->w1,mr_mip->w8);
            nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w1,mr_mip->w8);
            if (mr_mip->Asize!=0) nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w3,mr_mip->w8);
        }
    }
    else
    { /* a is not special */
        nres_modmult(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w3);
        nres_modmult(_MIPP_ mr_mip->w3,mr_mip->A,mr_mip->w3);
        nres_modmult(_MIPP_ p->X,p->X,mr_mip->w1);
        nres_modadd(_MIPP_ mr_mip->w1,mr_mip->w1,mr_mip->w8);
        nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w1,mr_mip->w8);
        nres_modadd(_MIPP_ mr_mip->w8,mr_mip->w3,mr_mip->w8);        
    }

/* w8 contains numerator of slope 3x^2+A.z^4  *
 * denominator is now placed in Z             */
    if (p->marker==MR_EPOINT_NORMALIZED)
        copy(p->Y,p->Z);
    else nres_modmult(_MIPP_ p->Y,p->Z,p->Z);

    nres_modadd(_MIPP_ p->Z,p->Z,p->Z);
    nres_modmult(_MIPP_ p->Y,p->Y,mr_mip->w2);
    nres_modmult(_MIPP_ p->X,mr_mip->w2,mr_mip->w3);
    nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w3,mr_mip->w3);
    nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w3,mr_mip->w3);
    nres_modmult(_MIPP_ mr_mip->w8,mr_mip->w8,mr_mip->w1);
    nres_modadd(_MIPP_ mr_mip->w3,mr_mip->w3,p->X);
    nres_modsub(_MIPP_ mr_mip->w1,p->X,p->X);
    nres_modmult(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w2);
    nres_modadd(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w2);
    nres_modadd(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w2);
    nres_modadd(_MIPP_ mr_mip->w2,mr_mip->w2,mr_mip->w2);
    nres_modsub(_MIPP_ mr_mip->w3,p->X,mr_mip->w3);
    nres_modmult(_MIPP_ mr_mip->w3,mr_mip->w8,mr_mip->w3);
    nres_modsub(_MIPP_ mr_mip->w3,mr_mip->w2,p->Y);
  
    p->marker=MR_EPOINT_GENERAL;
    return;
}
   
static BOOL ecurve_padd(_MIPD_ epoint *p,epoint *pa)
{ /* primitive add two epoints on the active ecurve - pa+=p;   *
   * note that if p is normalized, its Z coordinate isn't used */
 
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->coord==MR_AFFINE)
    {
        nres_modsub(_MIPP_ p->Y,pa->Y,mr_mip->w8);
        nres_modsub(_MIPP_ p->X,pa->X,mr_mip->w6);
        if (size(mr_mip->w6)==0)