mrecgf2m.c

大数运算库

    M=mr_mip->M;
    A=mr_mip->AA;
    B=mr_mip->BB;
    C=mr_mip->CC;
    convert(_MIPP_ 1,mr_mip->w1);     /* B */
    zero(mr_mip->w2);                 /* C */
    copy(x,mr_mip->w3);               /* F */
    copy(mr_mip->modulus,mr_mip->w4); /* G */
    forever
    {
        bits=zerobits(mr_mip->w3);
        shiftrightbits(mr_mip->w3,bits);
        shiftleftbits(mr_mip->w2,bits);
        k+=bits;    
        if (size(mr_mip->w3)==1) break;

        if (numbits(mr_mip->w3)<numbits(mr_mip->w4))
        { /* swap F & G, B & C */
            t=mr_mip->w3; mr_mip->w3=mr_mip->w4; mr_mip->w4=t;
            t=mr_mip->w1; mr_mip->w1=mr_mip->w2; mr_mip->w2=t;
        }
        add2(mr_mip->w3,mr_mip->w4,mr_mip->w3);
        add2(mr_mip->w1,mr_mip->w2,mr_mip->w1);
    }

    copy(mr_mip->w1,w);

    if (k==0) return TRUE;
    step=MIRACL;

    if (A<MIRACL) step=A;
    
    k1=1+M/MIRACL;       /* words from MSB to LSB */
    rs1=M%MIRACL;
    ls1=MIRACL-rs1;

    k2=1+A/MIRACL;   /* words from MSB to bit */
    rs2=A%MIRACL;
    ls2=MIRACL-rs2;

    if (B)
    { /* Pentanomial */
        if (C<MIRACL) step=C;

        k3=1+B/MIRACL;
        rs3=B%MIRACL;
        ls3=MIRACL-rs3;

        k4=1+C/MIRACL;
        rs4=C%MIRACL;
        ls4=MIRACL-rs4;
    }

    gw=w->w;
    while (k>0)
    {
        if (k>step) n=step;
        else        n=k;
 
        if (n==MIRACL) lsw=gw[0];
        else           lsw=gw[0]&(((mr_small)1<<n)-1);

        w->len=k1;
        if (rs1==0) gw[k1-1]^=lsw;
        else
        {
            w->len++;
            gw[k1]^=(lsw>>ls1);
            gw[k1-1]^=(lsw<<rs1);
        }
        if (rs2==0) gw[k2-1]^=lsw;
        else
        {
            gw[k2]^=(lsw>>ls2);
            gw[k2-1]^=(lsw<<rs2);
        }
        if (B)
        {
            if (rs3==0) gw[k3-1]^=lsw;
            else
            {
                gw[k3]^=(lsw>>ls3);
                gw[k3-1]^=(lsw<<rs3);
            }
            if (rs4==0) gw[k4-1]^=lsw;
            else
            {
                gw[k4]^=(lsw>>ls4);
                gw[k4-1]^=(lsw<<rs4);
            }
        }
        shiftrightbits(w,n);
        k-=n;
    }
    mr_lzero(w);
    return TRUE;
}

BOOL multi_inverse2(_MIPD_ int m,big *x,big *w)
{ /* find w[i]=1/x[i] mod f, for i=0 to m-1 */
    int i;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (m==0) return TRUE;
    if (m<0) return FALSE;

    if (x==w)
    {
        mr_berror(_MIPP_ MR_ERR_BAD_PARAMETERS);
        return FALSE;
    }
    if (m==1)
    {
        inverse2(_MIPP_ x[0],w[0]);
        return TRUE;
    }
    convert(_MIPP_ 1,w[0]);
    copy(x[0],w[1]);
    for (i=2;i<m;i++)
        modmult2(_MIPP_ w[i-1],x[i-1],w[i]);

    modmult2(_MIPP_ w[m-1],x[m-1],mr_mip->w6);
    if (size(mr_mip->w6)==0)
    {
        mr_berror(_MIPP_ MR_ERR_DIV_BY_ZERO);
        return FALSE;
    }

    inverse2(_MIPP_ mr_mip->w6,mr_mip->w6);  /* y=1/y */

    copy(x[m-1],mr_mip->w5);
    modmult2(_MIPP_ w[m-1],mr_mip->w6,w[m-1]);

    for (i=m-2;;i--)
    {
        if (i==0)
        {
            modmult2(_MIPP_ mr_mip->w5,mr_mip->w6,w[0]);
            break;
        }
        modmult2(_MIPP_ w[i],mr_mip->w5,w[i]);
        modmult2(_MIPP_ w[i],mr_mip->w6,w[i]);
        modmult2(_MIPP_ mr_mip->w5,x[i],mr_mip->w5);
    }
    return TRUE;
}

int trace2(_MIPD_ big x)
{
    int i;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    copy(x,mr_mip->w1);
    for (i=1;i<mr_mip->M;i++)
    {
        modsquare2(_MIPP_ mr_mip->w1,mr_mip->w1);
        add2(mr_mip->w1,x,mr_mip->w1);
    }    
    return (mr_mip->w1->w[0]&1);
}

void rand2(_MIPD_ big x)
{ /* random number */
    int i,k;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    zero(x);
    k=1+mr_mip->M/MIRACL;        
    x->len=k;
    for (i=0;i<k;i++) x->w[i]=brand(_MIPPO_ );
    mr_lzero(x);
    reduce2(_MIPP_ x,x);    
}

int parity2(big x)
{ /* return LSB */
   if (x->len==0) return 0;
   return (int)(x->w[0]%2);
}

BOOL quad2(_MIPD_ big b,big w)
{ /* Solves x^2 + x = b  for a root w  *
   * returns TRUE if a solution exists *
   * the "other" solution is w+1       */
    int i,M;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif

    M=mr_mip->M;
    copy(b,mr_mip->w1);
    if (M%2==1)
    { /* M is odd, so its the Half-Trace */
        copy(b,w);
        for (i=1;i<=(M-1)/2;i++)
        { 
            modsquare2(_MIPP_ w,w);
            modsquare2(_MIPP_ w,w);
            add2(w,mr_mip->w1,w);   
        } 
    }
    else
    {
        forever
        {
            rand2(_MIPP_ mr_mip->w2);
            zero(w);
            copy(mr_mip->w2,mr_mip->w3);
            for (i=1;i<M;i++)
            {
                modsquare2(_MIPP_ mr_mip->w3,mr_mip->w3);
                modmult2(_MIPP_ mr_mip->w3,mr_mip->w1,mr_mip->w4);
                modsquare2(_MIPP_ w,w);
                add2(w,mr_mip->w4,w);
                add2(mr_mip->w3,mr_mip->w2,mr_mip->w3);
            }    
            if (size(mr_mip->w3)!=0) break; 
        }
    }
    copy(w,mr_mip->w2);
    modsquare2(_MIPP_ mr_mip->w2,mr_mip->w2);
    add2(mr_mip->w2,w,mr_mip->w2);
    if (compare(mr_mip->w1,mr_mip->w2)==0) return TRUE;
    return FALSE;
}

void gf2m_dotprod(_MIPD_ int n,big *x,big *y,big w)
{ /* dot product - only one reduction! */
    int i;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    mr_mip->check=OFF;
    zero(mr_mip->w5);

    for (i=0;i<n;i++)
    {
        multiply2(_MIPP_ x[i],y[i],mr_mip->w0);
        add2(mr_mip->w5,mr_mip->w0,mr_mip->w5);
    }

    reduce2(_MIPP_ mr_mip->w5,mr_mip->w5);
    copy(mr_mip->w5,w);

    mr_mip->check=ON;
}

BOOL prepare_basis(_MIPD_ int m,int a,int b,int c,BOOL check)
{
    int i,k,sh;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return FALSE;

    if (b==0) c=0;
    if (m==mr_mip->M && a==mr_mip->AA && b==mr_mip->BB && c==mr_mip->CC)
        return TRUE;   /* its already prepared... */

    MR_IN(138)
    if (m <=0 || a<=0 || a>=m || b>=a) 
    {
        mr_berror(_MIPP_ MR_ERR_BAD_MODULUS);
        MR_OUT
        return FALSE;
    }
    
    mr_mip->M=m;
    mr_mip->AA=a;
    mr_mip->BB=0;
    mr_mip->CC=0;
    if (mr_mip->modulus==NULL) mr_mip->modulus=mirvar(_MIPP_ 0);
    else zero(mr_mip->modulus);

    k=1+m/MIRACL;
    mr_mip->modulus->len=k;
    sh=m%MIRACL;
    mr_mip->modulus->w[k-1]=((mr_small)1<<sh);
    mr_mip->modulus->w[0]^=1;
    mr_mip->modulus->w[a/MIRACL]^=((mr_small)1<<(a%MIRACL));
    if (b!=0)
    {
         mr_mip->BB=b;
         mr_mip->CC=c;
         mr_mip->modulus->w[b/MIRACL]^=((mr_small)1<<(b%MIRACL));
         mr_mip->modulus->w[c/MIRACL]^=((mr_small)1<<(c%MIRACL));
    }

    if (!check)
    {
        MR_OUT
        return TRUE;
    }

/* check for irreducibility of basis */

    zero(mr_mip->w4);
    mr_mip->w4->len=1;
    mr_mip->w4->w[0]=2;       /* f(t) = t */
    for (i=0;i<=m/2;i++)
    {
        modsquare2(_MIPP_ mr_mip->w4,mr_mip->w4);
        incr2(mr_mip->w4,2,mr_mip->w5);
        gcd2(_MIPP_ mr_mip->w5,mr_mip->modulus,mr_mip->w5);
        if (size(mr_mip->w5)!=1)
        {
            mr_berror(_MIPP_ MR_ERR_NOT_IRREDUC);
            MR_OUT
            return FALSE;
        }
    }
                   
    MR_OUT
    return TRUE;
}

/* Initialise with Trinomial or Pentanomial      *
 * t^m  + t^a + 1 OR t^m + t^a +t^b + t^c + 1    *
 * Set b=0 for pentanomial. a2 is usually 0 or 1 *
 * m negative indicates a super-singular curve   */

BOOL ecurve2_init(_MIPD_ int m,int a,int b,int c,big a2,big a6,BOOL check,int type)
{
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    int i;

/* catch some nonsense conditions */

    if (mr_mip->ERNUM) return FALSE;

    mr_mip->SS=FALSE;
    if (m<0)
    { /* its a supersingular curve! */
        mr_mip->SS=TRUE;
        type=MR_AFFINE;     /* always AFFINE */
        m=-m;
        if (size(a2)!=1) return FALSE;
        if (size(a6) >1) return FALSE;
    }
    if (size(a2)<0) return FALSE;
    if (size(a6)<0) return FALSE;

    MR_IN(123)

    if (!prepare_basis(_MIPP_ m,a,b,c,check))
    { /* unable to set the basis */
        MR_OUT
        return FALSE;
    }    

    mr_mip->Asize=size(a2);    
    mr_mip->Bsize=size(a6);

    if (mr_mip->Asize==MR_TOOBIG)
    { 
        if (mr_mip->A==NULL) mr_mip->A=mirvar(_MIPP_ 0);
        copy(a2,mr_mip->A);
    }

    if (mr_mip->Bsize==MR_TOOBIG)
    { 
        if (mr_mip->B==NULL) mr_mip->B=mirvar(_MIPP_ 0);
        copy(a6,mr_mip->B);
    }

 /* Use C to store B^(2^M-2) - required for projective doubling */
    if (!mr_mip->SS)
    {
        if (mr_mip->C==NULL) mr_mip->C=mirvar(_MIPP_ 0);
        copy(a6,mr_mip->C);
        for (i=1;i<m-1;i++) modsquare2(_MIPP_ mr_mip->C,mr_mip->C);
    }

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

epoint* epoint2_init(_MIPDO_ )
{
    epoint *p;
    char *ptr;
#ifndef MR_GENERIC_MT
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return NULL;
    
    MR_IN(124)

/* 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 epoint2_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 epoint2_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/x     *
   * (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(125)

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

    valid=FALSE;       

    if (mr_mip->SS)
    { /* Super-singular - calculate x^3+x+B */
        copy (x,p->X);
        modsquare2(_MIPP_ p->X,mr_mip->w5);           /* w5=x^2 */
        modmult2(_MIPP_ mr_mip->w5,p->X,mr_mip->w5);  /* w5=x^3 */
        add2(mr_mip->w5,p->X,mr_mip->w5);             
        incr2(mr_mip->w5,mr_mip->Bsize,mr_mip->w5);  /* w5=x^3+x+B */
        if (x!=y)
        { /* compare with y^2+y */
            copy(y,p->Y);
            modsquare2(_MIPP_ p->Y,mr_mip->w1);
            add2(mr_mip->w1,p->Y,mr_mip->w1);
            if (compare(mr_mip->w1,mr_mip->w5)==0) valid=TRUE;
        }
        else
        { /* no y supplied - calculate one. Solve quadratic */
            valid=quad2(_MIPP_ mr_mip->w5,mr_mip->w5);
            incr2(mr_mip->w5,cb^parity2(mr_mip->w5),p->Y);
        }
    } 
    else
    { /* calculate x^3+Ax^2+B */
        copy(x,p->X);

        modsquare2(_MIPP_ p->X,mr_mip->w6);           /* w6=x^2 */
        modmult2(_MIPP_ mr_mip->w6,p->X,mr_mip->w5);  /* w5=x^3 */

        if (mr_mip->Asize==MR_TOOBIG)
            copy(mr_mip->A,mr_mip->w1);
        else
            convert(_MIPP_ mr_mip->Asize,mr_mip->w1);
        modmult2(_MIPP_ mr_mip->w6,mr_mip->w1,mr_mip->w0);
        add2(mr_mip->w5,mr_mip->w0,mr_mip->w5);

        if (mr_mip->Bsize==MR_TOOBIG)
            add2(mr_mip->w5,mr_mip->B,mr_mip->w5);    /* w5=x^3+Ax^2+B */
        else
            incr2(mr_mip->w5,mr_mip->Bsize,mr_mip->w5); 
        if (x!=y)
        { /* compare with y^2+xy */
            copy(y,p->Y);
            modsquare2(_MIPP_ p->Y,mr_mip->w2);