mrecn2_opt.c

miracl大数库 miracl大数库 miracl大数库

#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

#ifndef MR_STATIC
    char *mem = memalloc(_MIPP_ MR_DOS_2);
#else
    char mem[MR_BIG_RESERVE(MR_DOS_2)];       
 	memset(mem, 0, MR_BIG_RESERVE(MR_DOS_2));   
#endif

    twist=mr_mip->TWIST;
    j=0;
    sz=calc_n(win);

    A.a= mirvar_mem(_MIPP_ mem, j++);
    A.b= mirvar_mem(_MIPP_ mem, j++);
    B.a= mirvar_mem(_MIPP_ mem, j++);
    B.b= mirvar_mem(_MIPP_ mem, j++);
    C.a= mirvar_mem(_MIPP_ mem, j++);
    C.b= mirvar_mem(_MIPP_ mem, j++);
    D.a= mirvar_mem(_MIPP_ mem, j++);
    D.b= mirvar_mem(_MIPP_ mem, j++);
    E.a= mirvar_mem(_MIPP_ mem, j++);
    E.b= mirvar_mem(_MIPP_ mem, j++);
    T.a= mirvar_mem(_MIPP_ mem, j++);
    T.b= mirvar_mem(_MIPP_ mem, j++);
    W.a= mirvar_mem(_MIPP_ mem, j++);
    W.b= mirvar_mem(_MIPP_ mem, j++);

    for (i=0;i<sz;i++)
    {
        d[i].a= mirvar_mem(_MIPP_ mem, j++);
        d[i].b= mirvar_mem(_MIPP_ mem, j++);
        e[i].a= mirvar_mem(_MIPP_ mem, j++);
        e[i].b= mirvar_mem(_MIPP_ mem, j++);
    }

    zzn2_add_i(&(PT[0].y),&(PT[0].y),&d[0]);   /* 1. d_0=2.y */
    zzn2_sqr_i(&d[0],&C);                      /* 2. C=d_0^2 */

    zzn2_sqr_i(&(PT[0].x),&T);
    zzn2_add_i(&T,&T,&A);
    zzn2_add_i(&T,&A,&T);
           
    if (mr_abs(mr_mip->Asize)<MR_TOOBIG) zzn2_from_int(_MIPP_ mr_mip->Asize,&A);
    else zzn2_from_zzn(mr_mip->A,&A);
        
    if (twist)
    {
        zzn2_txx_i(&A);
        zzn2_txx_i(&A);
    }
    zzn2_add_i(&A,&T,&A);             /* 3. A=3x^2+a */
    zzn2_copy_i(&A,&W);

    zzn2_add_i(&C,&C,&B);
    zzn2_add_i(&B,&C,&B);
    zzn2_mul_i(&B,&(PT[0].x),&B);     /* 4. B=3C.x */

    zzn2_sqr_i(&A,&d[1]);
    zzn2_sub_i(&d[1],&B,&d[1]);       /* 5. d_1=A^2-B */

    zzn2_sqr_i(&d[1],&E);             /* 6. E=d_1^2 */
    
    zzn2_mul_i(&B,&E,&B);             /* 7. B=E.B */

    zzn2_sqr_i(&C,&C);                /* 8. C=C^2 */

    zzn2_mul_i(&E,&d[1],&D);          /* 9. D=E.d_1 */

    zzn2_mul_i(&A,&d[1],&A);
    zzn2_add_i(&A,&C,&A);
    zzn2_negate(_MIPP_ &A,&A);             /* 10. A=-d_1*A-C */

    zzn2_add_i(&D,&D,&T);
    zzn2_sqr_i(&A,&d[2]);
    zzn2_sub_i(&d[2],&T,&d[2]);
    zzn2_sub_i(&d[2],&B,&d[2]);       /* 11. d_2=A^2-2D-B */

    if (sz>3)
    {
        zzn2_sqr_i(&d[2],&E);             /* 12. E=d_2^2 */

        zzn2_add_i(&T,&D,&T);
        zzn2_add_i(&T,&B,&T);
        zzn2_mul_i(&T,&E,&B);             /* 13. B=E(B+3D) */
        
        zzn2_add_i(&A,&A,&T);
        zzn2_add_i(&C,&T,&C);
        zzn2_mul_i(&C,&D,&C);             /* 14. C=D(2A+C) */

        zzn2_mul_i(&d[2],&E,&D);          /* 15. D=E.d_2 */

        zzn2_mul_i(&A,&d[2],&A);
        zzn2_add_i(&A,&C,&A);
        zzn2_negate(_MIPP_ &A,&A);             /* 16. A=-d_2*A-C */

 
        zzn2_sqr_i(&A,&d[3]);
        zzn2_sub_i(&d[3],&D,&d[3]);
        zzn2_sub_i(&d[3],&B,&d[3]);       /* 17. d_3=A^2-D-B */

        for (i=4;i<sz;i++)
        {
            zzn2_sqr_i(&d[i-1],&E);       /* 19. E=d(i-1)^2 */
            zzn2_mul_i(&B,&E,&B);         /* 20. B=E.B */
            zzn2_mul_i(&C,&D,&C);         /* 21. C=D.C */
            zzn2_mul_i(&E,&d[i-1],&D);    /* 22. D=E.d(i-1) */

            zzn2_mul_i(&A,&d[i-1],&A);
            zzn2_add_i(&A,&C,&A);
            zzn2_negate(_MIPP_ &A,&A);         /* 23. A=-d(i-1)*A-C */

            zzn2_sqr_i(&A,&d[i]);
            zzn2_sub_i(&d[i],&D,&d[i]);
            zzn2_sub_i(&d[i],&B,&d[i]);   /* 24. d(i)=A^2-D-B */
        }
    }

    zzn2_copy_i(&d[0],&e[0]);
    for (i=1;i<sz;i++)
        zzn2_mul_i(&e[i-1],&d[i],&e[i]);
       
    zzn2_copy_i(&e[sz-1],&A);
    zzn2_inv_i(_MIPP_ &A);

    for (i=sz-1;i>0;i--)
    {
        zzn2_copy_i(&d[i],&B);
        zzn2_mul_i(&e[i-1],&A,&d[i]);  
        zzn2_mul_i(&A,&B,&A);
    }
    zzn2_copy_i(&A,&d[0]);

    for (i=1;i<sz;i++)
    {
        zzn2_sqr_i(&e[i-1],&T);
        zzn2_mul_i(&d[i],&T,&d[i]); /** */
    }

    zzn2_mul_i(&W,&d[0],&W);
    zzn2_sqr_i(&W,&A);
    zzn2_sub_i(&A,&(PT[0].x),&A);
    zzn2_sub_i(&A,&(PT[0].x),&A);
    zzn2_sub_i(&(PT[0].x),&A,&B);
    zzn2_mul_i(&B,&W,&B);
    zzn2_sub_i(&B,&(PT[0].y),&B);

    zzn2_sub_i(&B,&(PT[0].y),&T);
    zzn2_mul_i(&T,&d[1],&T);

    zzn2_sqr_i(&T,&(PT[1].x));
    zzn2_sub_i(&(PT[1].x),&A,&(PT[1].x));
    zzn2_sub_i(&(PT[1].x),&(PT[0].x),&(PT[1].x));

    zzn2_sub_i(&A,&(PT[1].x),&(PT[1].y));
    zzn2_mul_i(&(PT[1].y),&T,&(PT[1].y));
    zzn2_sub_i(&(PT[1].y),&B,&(PT[1].y));

    for (i=2;i<sz;i++)
    {
        zzn2_sub_i(&(PT[i-1].y),&B,&T);
        zzn2_mul_i(&T,&d[i],&T);

        zzn2_sqr_i(&T,&(PT[i].x));
        zzn2_sub_i(&(PT[i].x),&A,&(PT[i].x));
        zzn2_sub_i(&(PT[i].x),&(PT[i-1].x),&(PT[i].x));

        zzn2_sub_i(&A,&(PT[i].x),&(PT[i].y));
        zzn2_mul_i(&(PT[i].y),&T,&(PT[i].y));
        zzn2_sub_i(&(PT[i].y),&B,&(PT[i].y));
    }
    for (i=0;i<sz;i++) PT[i].marker=MR_EPOINT_NORMALIZED;

#ifndef MR_STATIC
    memkill(_MIPP_ mem, MR_DOS_2);
#else
    memset(mem, 0, MR_BIG_RESERVE(MR_DOS_2));
#endif
}

#ifndef MR_DOUBLE_BIG
#define MR_MUL_RESERVE (1+4*MR_STR_SZ_2)
#else
#define MR_MUL_RESERVE (2+4*MR_STR_SZ_2)
#endif

int ecn2_mul(_MIPD_ big k,ecn2 *P)
{
    int i,j,nb,n,nbs,nzs,nadds;
    big h;
    ecn2 T[MR_STR_SZ_2];
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

#ifndef MR_STATIC
    char *mem = memalloc(_MIPP_ MR_MUL_RESERVE);
#else
    char mem[MR_BIG_RESERVE(MR_MUL_RESERVE)];
    memset(mem, 0, MR_BIG_RESERVE(MR_MUL_RESERVE));
#endif

    j=0;
#ifndef MR_DOUBLE_BIG
    h=mirvar_mem(_MIPP_ mem, j++);
#else
    h=mirvar_mem(_MIPP_ mem, j); j+=2;
#endif
    for (i=0;i<MR_STR_SZ_2;i++)
    {
        T[i].x.a= mirvar_mem(_MIPP_ mem, j++);
        T[i].x.b= mirvar_mem(_MIPP_ mem, j++);
        T[i].y.a= mirvar_mem(_MIPP_ mem, j++);
        T[i].y.b= mirvar_mem(_MIPP_ mem, j++);
    }

    MR_IN(207)

    ecn2_norm(_MIPP_ P);

	nadds=0;
    premult(_MIPP_ k,3,h);
    ecn2_copy(P,&T[0]);
    ecn2_dos(_MIPP_ MR_WIN_SZ_2,T);
    nb=logb2(_MIPP_ h);

    for (i=nb-2;i>=1;)
    {
        if (mr_mip->user!=NULL) (*mr_mip->user)();
        n=mr_naf_window(_MIPP_ k,h,i,&nbs,&nzs,MR_WIN_SZ_2);
 
        for (j=0;j<nbs;j++) ecn2_add(_MIPP_ P,P);
       
        if (n>0) {nadds++; ecn2_add(_MIPP_ &T[n/2],P);}
        if (n<0) {nadds++; ecn2_sub(_MIPP_ &T[(-n)/2],P);}
        i-=nbs;
        if (nzs)
        {
            for (j=0;j<nzs;j++) ecn2_add(_MIPP_ P,P);
            i-=nzs;
        }
    }

    ecn2_norm(_MIPP_ P);

    MR_OUT

#ifndef MR_STATIC
    memkill(_MIPP_ mem, MR_MUL_RESERVE);
#else
    memset(mem, 0, MR_BIG_RESERVE(MR_MUL_RESERVE));
#endif
	return nadds;
}

/* Double addition, using Joint Sparse Form */
/* R=aP+bQ */

#define MR_MUL2_JSF_RESERVE 20

int ecn2_mul2_jsf(_MIPD_ big a,ecn2 *P,big b,ecn2 *Q,ecn2 *R)
{
    int e1,h1,e2,h2,bb,nadds;
    ecn2 P1,P2,PS,PD;
    big c,d,e,f;
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

#ifndef MR_STATIC
    char *mem = memalloc(_MIPP_ MR_MUL2_JSF_RESERVE);
#else
    char mem[MR_BIG_RESERVE(MR_MUL2_JSF_RESERVE)];
    memset(mem, 0, MR_BIG_RESERVE(MR_MUL2_JSF_RESERVE));
#endif

    c = mirvar_mem(_MIPP_ mem, 0);
    d = mirvar_mem(_MIPP_ mem, 1);
    e = mirvar_mem(_MIPP_ mem, 2);
    f = mirvar_mem(_MIPP_ mem, 3);
    P1.x.a= mirvar_mem(_MIPP_ mem, 4);
    P1.x.b= mirvar_mem(_MIPP_ mem, 5);
    P1.y.a= mirvar_mem(_MIPP_ mem, 6);
    P1.y.b= mirvar_mem(_MIPP_ mem, 7);
    P2.x.a= mirvar_mem(_MIPP_ mem, 8);
    P2.x.b= mirvar_mem(_MIPP_ mem, 9);
    P2.y.a= mirvar_mem(_MIPP_ mem, 10);
    P2.y.b= mirvar_mem(_MIPP_ mem, 11);
    PS.x.a= mirvar_mem(_MIPP_ mem, 12);
    PS.x.b= mirvar_mem(_MIPP_ mem, 13);
    PS.y.a= mirvar_mem(_MIPP_ mem, 14);
    PS.y.b= mirvar_mem(_MIPP_ mem, 15);
    PD.x.a= mirvar_mem(_MIPP_ mem, 16);
    PD.x.b= mirvar_mem(_MIPP_ mem, 17);
    PD.y.a= mirvar_mem(_MIPP_ mem, 18);
    PD.y.b= mirvar_mem(_MIPP_ mem, 19);

    MR_IN(206)

    ecn2_norm(_MIPP_ Q); 
    ecn2_copy(Q,&P2); 

    copy(b,d);
    if (size(d)<0) 
    {
        negify(d,d);
        ecn2_negate(_MIPP_ &P2,&P2);
    }

    ecn2_norm(_MIPP_ P); 
    ecn2_copy(P,&P1); 

    copy(a,c);
    if (size(c)<0) 
    {
        negify(c,c);
        ecn2_negate(_MIPP_ &P1,&P1);
    }

    mr_jsf(_MIPP_ d,c,e,d,f,c);   /* calculate joint sparse form */
 
    if (compare(e,f)>0) bb=logb2(_MIPP_ e)-1;
    else                bb=logb2(_MIPP_ f)-1;

    ecn2_add_sub(_MIPP_ &P1,&P2,&PS,&PD);
    ecn2_zero(R);
	nadds=0;
    while (bb>=0) 
    { /* add/subtract method */
        if (mr_mip->user!=NULL) (*mr_mip->user)();
        ecn2_add(_MIPP_ R,R);
        e1=h1=e2=h2=0;

        if (mr_testbit(_MIPP_ d,bb)) e2=1;
        if (mr_testbit(_MIPP_ e,bb)) h2=1;
        if (mr_testbit(_MIPP_ c,bb)) e1=1;
        if (mr_testbit(_MIPP_ f,bb)) h1=1;

        if (e1!=h1)
        {
            if (e2==h2)
            {
                if (h1==1) {ecn2_add(_MIPP_ &P1,R); nadds++;}
                else       {ecn2_sub(_MIPP_ &P1,R); nadds++;}
            }
            else
            {
                if (h1==1)
                {
                    if (h2==1) {ecn2_add(_MIPP_ &PS,R); nadds++;}
                    else       {ecn2_add(_MIPP_ &PD,R); nadds++;}
                }
                else
                {
                    if (h2==1) {ecn2_sub(_MIPP_ &PD,R); nadds++;}
                    else       {ecn2_sub(_MIPP_ &PS,R); nadds++;}
                }
            }
        }
        else if (e2!=h2)
        {
            if (h2==1) {ecn2_add(_MIPP_ &P2,R); nadds++;}
            else       {ecn2_sub(_MIPP_ &P2,R); nadds++;}
        }
        bb-=1;
    }
    ecn2_norm(_MIPP_ R); 

    MR_OUT
#ifndef MR_STATIC
    memkill(_MIPP_ mem, MR_MUL2_JSF_RESERVE);
#else
    memset(mem, 0, MR_BIG_RESERVE(MR_MUL2_JSF_RESERVE));
#endif
	return nadds;

}

/* General purpose multi-exponentiation engine, using inter-leaving algorithm. Calculate aP+bQ+cR+dS...
   Inputs are divided into two groups of sizes wa<4 and wb<4. For the first group if the points are fixed the 
   first precomputed Table Ta[] may be taken from ROM. For the second group if the points are variable Tb[j] will
   have to computed online. Each group has its own window size, wina (=5?) and winb (=4?) respectively. The values 
   a,b,c.. are provided in ma[] and mb[], and 3.a,3.b,3.c (as required by the NAF) are provided in ma3[] and 
   mb3[]. If only one group is required, set wb=0 and pass NULL pointers.
   */

int ecn2_muln_engine(_MIPD_ int wa,int wina,int wb,int winb,big *ma,big *ma3,big *mb,big *mb3,ecn2 *Ta,ecn2 *Tb,ecn2 *R)
{ /* general purpose interleaving algorithm engine for multi-exp */
    int i,j,tba[4],pba[4],na[4],sa[4],tbb[4],pbb[4],nb[4],sb[4],nbits,nbs,nzs;
    int sza,szb,nadds;
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

    sza=calc_n(wina);
    szb=calc_n(winb);

    ecn2_zero(R);

    nbits=0;
    for (i=0;i<wa;i++) {sa[i]=exsign(ma[i]); tba[i]=0; j=logb2(_MIPP_ ma3[i]); if (j>nbits) nbits=j; }
    for (i=0;i<wb;i++) {sb[i]=exsign(mb[i]); tbb[i]=0; j=logb2(_MIPP_ mb3[i]); if (j>nbits) nbits=j; }
    
    nadds=0;
    for (i=nbits-1;i>=1;i--)
    {
        if (mr_mip->user!=NULL) (*mr_mip->user)();
        if (R->marker!=MR_EPOINT_INFINITY) ecn2_add(_MIPP_ R,R);
        for (j=0;j<wa;j++)
        { /* deal with the first group */
            if (tba[j]==0)
            {
                na[j]=mr_naf_window(_MIPP_ ma[j],ma3[j],i,&nbs,&nzs,wina);
                tba[j]=nbs+nzs;
                pba[j]=nbs;
            }
            tba[j]--;  pba[j]--; 
            if (pba[j]==0)
            {
                if (sa[j]==PLUS)
                {
                    if (na[j]>0) {ecn2_add(_MIPP_ &Ta[j*sza+na[j]/2],R); nadds++;}
                    if (na[j]<0) {ecn2_sub(_MIPP_ &Ta[j*sza+(-na[j])/2],R); nadds++;}
                }
                else
                {
                    if (na[j]>0) {ecn2_sub(_MIPP_ &Ta[j*sza+na[j]/2],R); nadds++;}
                    if (na[j]<0) {ecn2_add(_MIPP_ &Ta[j*sza+(-na[j])/2],R); nadds++;}
                }
            }         
        }
        for (j=0;j<wb;j++)
        { /* deal with the second group */
            if (tbb[j]==0)
            {
                nb[j]=mr_naf_window(_MIPP_ mb[j],mb3[j],i,&nbs,&nzs,winb);
                tbb[j]=nbs+nzs;
                pbb[j]=nbs;
            }
            tbb[j]--;  pbb[j]--; 
            if (pbb[j]==0)
            {
                if (sb[j]==PLUS)
                {
                    if (nb[j]>0) {ecn2_add(_MIPP_ &Tb[j*szb+nb[j]/2],R);  nadds++;}
                    if (nb[j]<0) {ecn2_sub(_MIPP_ &Tb[j*szb+(-nb[j])/2],R);  nadds++;}
                }
                else
                {
                    if (nb[j]>0) {ecn2_sub(_MIPP_ &Tb[j*szb+nb[j]/2],R);  nadds++;}
                    if (nb[j]<0) {ecn2_add(_MIPP_ &Tb[j*szb+(-nb[j])/2],R);  nadds++;}
                }
            }         
        }
    }
    ecn2_norm(_MIPP_ R);  
    return nadds;
}

/* Routines to support Galbraith, Lin, Scott (GLS) method for ECC */
/* requires an endomorphism psi */

/* *********************** */

/* Precompute T - first half from i.P, second half from i.psi(P) */ 

void ecn2_precomp_gls(_MIPD_ int win,ecn2 *P,zzn2 *psi,ecn2 *T)
{
    int i,j,sz;

#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

    j=0;
    sz=calc_n(win);

    MR_IN(219)

    ecn2_norm(_MIPP_ P);
    ecn2_copy(P,&T[0]);
    
    ecn2_dos(_MIPP_ win,T); /* precompute table */

    for (i=sz;i<sz+sz;i++)
    {
        ecn2_copy(&T[i-sz],&T[i]);
        ecn2_psi(_MIPP_ psi,&T[i]);
    }

    MR_OUT
}

/* Calculate a[0].P+a[1].psi(P) using interleaving method */

#define MR_MUL2_GLS_RESERVE (2+2*MR_STR_SZ_2*4)

int ecn2_mul2_gls(_MIPD_ big *a,ecn2 *P,zzn2 *psi,ecn2 *R)
{
    int i,j,nadds;
    ecn2 T[2*MR_STR_SZ_2];
    big a3[2];
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

#ifndef MR_STATIC
    char *mem = memalloc(_MIPP_ MR_MUL2_GLS_RESERVE);
#else