mrarth1.c

比较新的功能强大的rsa算法源代码,方便使用.

/*
 *   MIRACL arithmetic routines 1 - multiplying and dividing BIG NUMBERS by  
 *   integer numbers.
 *   mrarth1.c
 *
 *   Copyright (c) 1988-1997 Shamus Software Ltd.
 */

#include "miracl.h"

#ifdef MR_FP
#include <math.h>
#endif

#ifdef MR_FP_ROUNDING
#ifdef __GNUC__
#include <ieeefp.h>
#endif

/* Invert n and set FP rounding. 
 * Set to round up
 * Calculate 1/n
 * set to round down (towards zero)
 * If rounding cannot be controlled, this function returns 0.0 */

mr_large mr_invert(mr_small n)
{
    mr_large inn;
    int up=  0x1BFF; 

#ifdef _MSC_VER
  #ifdef MR_NOASM
#define NO_EXTENDED
  #endif 
#endif

#ifdef NO_EXTENDED
    int down=0x1EFF;
#else
    int down=0x1FFF;
#endif

#ifdef __TURBOC__
    asm
    {
        fldcw WORD PTR up
        fld1
        fld QWORD PTR n;
        fdiv
        fstp TBYTE PTR inn;
        fldcw WORD PTR down;
    }
    return inn;   
#endif
#ifdef _MSC_VER
    _asm
    {
        fldcw WORD PTR up
        fld1
        fld QWORD PTR n;
        fdiv
        fstp QWORD  PTR inn;
        fldcw WORD PTR down;
    }
    return inn;   
#endif
#ifdef __GNUC__
#ifdef i386
    __asm__ __volatile__ (
    "fldcw %2\n"
    "fld1\n"
    "fldl %1\n"
    "fdivrp\n"
    "fstpt %0\n"
    "fldcw %3\n"
    : "=m"(inn)
    : "m"(n),"m"(up),"m"(down)
    : "memory"
    );
    return inn;   
#else
    fpsetround(FP_RP);
    inn=(mr_large)1.0/n;
    fpsetround(FP_RZ);
    return inn;
#endif
#endif
    return 0.0L;   
}

#endif

void mr_pmul(_MIPD_ big x,mr_small sn,big z)
{ 
    int m,xl;
    mr_lentype sx;
    mr_small carry,*xg,*zg;

#ifdef MR_ITANIUM
    mr_small tm;
#endif
#ifdef MR_WIN64
    mr_small tm;
#endif
#ifdef MR_NOASM
    union doubleword dble;
    mr_large dbled;
    mr_large ldres;
#endif
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif

    if (x!=z)
    {
        zero(z);
        if (sn==0) return;
    }
    else if (sn==0)
    {
        zero(z);
        return;
    }
    m=0;
    carry=0;
    sx=x->len&MR_MSBIT;
    xl=(int)(x->len&MR_OBITS);

#ifndef MR_SIMPLE_BASE
    if (mr_mip->base==0) 
    {
#endif
#ifndef MR_NOFULLWIDTH
        xg=x->w; zg=z->w;
/* inline 8086 assembly - substitutes for loop below */
#ifdef INLINE_ASM
#if INLINE_ASM == 1
        ASM cld
        ASM mov cx,xl
        ASM or cx,cx
        ASM je out1
#ifdef MR_LMM
        ASM push ds
        ASM push es
        ASM les di,DWORD PTR zg
        ASM lds si,DWORD PTR xg
#else
        ASM mov ax,ds
        ASM mov es,ax
        ASM mov di,zg
        ASM mov si,xg
#endif
        ASM mov bx,sn
        ASM push bp
        ASM xor bp,bp
    tcl1:
        ASM lodsw
        ASM mul bx
        ASM add ax,bp
        ASM adc dx,0
        ASM stosw
        ASM mov bp,dx
        ASM loop tcl1

        ASM mov ax,bp
        ASM pop bp
#ifdef MR_LMM
        ASM pop es
        ASM pop ds
#endif
        ASM mov carry,ax
     out1: 
#endif
#if INLINE_ASM == 2
        ASM cld
        ASM mov cx,xl
        ASM or cx,cx
        ASM je out1
#ifdef MR_LMM
        ASM push ds
        ASM push es
        ASM les di,DWORD PTR zg
        ASM lds si,DWORD PTR xg
#else
        ASM mov ax,ds
        ASM mov es,ax
        ASM mov di,zg
        ASM mov si,xg
#endif
        ASM mov ebx,sn
        ASM push ebp
        ASM xor ebp,ebp
    tcl1:
        ASM lodsd
        ASM mul ebx
        ASM add eax,ebp
        ASM adc edx,0
        ASM stosd
        ASM mov ebp,edx
        ASM loop tcl1

        ASM mov eax,ebp
        ASM pop ebp
#ifdef MR_LMM
        ASM pop es
        ASM pop ds
#endif
        ASM mov carry,eax
     out1: 
#endif
#if INLINE_ASM == 3
        ASM mov ecx,xl
        ASM or ecx,ecx
        ASM je out1
        ASM mov ebx,sn
        ASM mov edi,zg
        ASM mov esi,xg
        ASM push ebp
        ASM xor ebp,ebp
    tcl1:
        ASM mov eax,[esi]
        ASM add esi,4
        ASM mul ebx
        ASM add eax,ebp
        ASM adc edx,0
        ASM mov [edi],eax
        ASM add edi,4
        ASM mov ebp,edx
        ASM dec ecx
        ASM jnz tcl1

        ASM mov eax,ebp
        ASM pop ebp
        ASM mov carry,eax
     out1: 
#endif
#if INLINE_ASM == 4

        ASM (
           "movl %4,%%ecx\n"
           "orl  %%ecx,%%ecx\n"
           "je 1f\n"
           "movl %3,%%ebx\n"
           "movl %1,%%edi\n"
           "movl %2,%%esi\n"
           "pushl %%ebp\n"
           "xorl %%ebp,%%ebp\n"  
        "0:\n"  
           "movl (%%esi),%%eax\n"
           "addl $4,%%esi\n"
           "mull %%ebx\n"
           "addl %%ebp,%%eax\n"
           "adcl $0,%%edx\n"
           "movl %%eax,(%%edi)\n"
           "addl $4,%%edi\n"
           "movl %%edx,%%ebp\n"
           "decl %%ecx\n"
           "jnz 0b\n"
 
           "movl %%ebp,%%eax\n"
           "popl %%ebp\n"
           "movl %%eax,%0\n"
        "1:"  
        :"=m"(carry)
        :"m"(zg),"m"(xg),"m"(sn),"m"(xl)
        :"eax","edi","esi","ebx","ecx","edx","memory"
        );

#endif
#endif
#ifndef INLINE_ASM
        for (m=0;m<xl;m++)
#ifdef MR_NOASM
        {
            dble.d=(mr_large)x->w[m]*sn+carry;
            carry=dble.h[MR_TOP];
            z->w[m]=dble.h[MR_BOT];
        }
#else
            carry=muldvd(x->w[m],sn,carry,&z->w[m]);
#endif
#endif
        if (carry>0)
        {
            m=xl;
            if (m>=mr_mip->nib && mr_mip->check)
            {
                mr_berror(_MIPP_ MR_ERR_OVERFLOW);
                return;
            }
            z->w[m]=carry;
            z->len=m+1;
        }
        else z->len=xl;
#endif
#ifndef MR_SIMPLE_BASE
    }
    else while (m<xl || carry>0)
    { /* multiply each digit of x by n */ 
    
        if (m>mr_mip->nib && mr_mip->check)
        {
            mr_berror(_MIPP_ MR_ERR_OVERFLOW);
            return;
        }
#ifdef MR_NOASM
        dbled=(mr_large)x->w[m]*sn+carry;
 #ifdef MR_FP_ROUNDING
        carry=(mr_small)MR_LROUND(dbled*mr_mip->inverse_base);
 #else
  #ifndef MR_FP
        if (mr_mip->base==mr_mip->base2)
          carry=(mr_small)(dbled>>mr_mip->lg2b);
        else 
  #endif  
          carry=(mr_small)MR_LROUND(dbled/mr_mip->base);
 #endif
        z->w[m]=(mr_small)(dbled-(mr_large)carry*mr_mip->base);
#else
 #ifdef MR_FP_ROUNDING
        carry=imuldiv(x->w[m],sn,carry,mr_mip->base,mr_mip->inverse_base,&z->w[m]);
 #else
        carry=muldiv(x->w[m],sn,carry,mr_mip->base,&z->w[m]);
 #endif
#endif

        m++;
        z->len=m;
    }
#endif
    if (z->len!=0) z->len|=sx;
}

void premult(_MIPD_ big x,int n,big z)
{ /* premultiply a big number by an int z=x.n */
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif
    if (mr_mip->ERNUM) return;

    MR_IN(9)


#ifdef MR_FLASH
    if (mr_notint(x))
    {
        mr_berror(_MIPP_ MR_ERR_INT_OP);
        MR_OUT
        return;
    }
#endif
    if (n==0)  /* test for some special cases  */
    {
        zero(z);
        MR_OUT
        return;
    }
    if (n==1)
    {
        copy(x,z);
        MR_OUT
        return;
    }
    if (n<0)
    {
        n=(-n);
        mr_pmul(_MIPP_ x,(mr_small)n,z);
        if (z->len!=0) z->len^=MR_MSBIT;
    }
    else mr_pmul(_MIPP_ x,(mr_small)n,z);
    MR_OUT
}

#ifdef MR_FP_ROUNDING
mr_small mr_sdiv(_MIPD_ big x,mr_small sn,mr_large isn,big z)
#else
mr_small mr_sdiv(_MIPD_ big x,mr_small sn,big z)
#endif
{
    int i,xl;
    mr_small sr,*xg,*zg;
#ifdef MR_NOASM
    union doubleword dble;
    mr_large dbled;
    mr_large ldres;
#endif
#ifdef MR_OS_THREADS
    miracl *mr_mip=get_mip();
#endif
    sr=0;
    xl=(int)(x->len&MR_OBITS);
    if (x!=z) zero(z);
#ifndef MR_SIMPLE_BASE
    if (mr_mip->base==0) 
    {
#endif
#ifndef MR_NOFULLWIDTH
        xg=x->w; zg=z->w;
/* inline - substitutes for loop below */
#ifdef INLINE_ASM
#if INLINE_ASM == 1
        ASM std
        ASM mov cx,xl
        ASM or cx,cx
        ASM je out2
        ASM mov bx,cx
        ASM shl bx,1
        ASM sub bx,2
#ifdef MR_LMM
        ASM push ds
        ASM push es
        ASM les di,DWORD PTR zg
        ASM lds si,DWORD PTR xg
#else
        ASM mov ax,ds
        ASM mov es,ax
        ASM mov di,zg
        ASM mov si,xg
#endif
        ASM add si,bx
        ASM add di,bx
        ASM mov bx,sn
        ASM push bp
        ASM xor bp,bp
    tcl2:
        ASM mov dx,bp
        ASM lodsw
        ASM div bx
        ASM mov bp,dx
        ASM stosw
        ASM loop tcl2

        ASM mov ax,bp
        ASM pop bp
#ifdef MR_LMM
        ASM pop es
        ASM pop ds
#endif
        ASM mov sr,ax
     out2:
        ASM cld
#endif
#if INLINE_ASM == 2
        ASM std
        ASM mov cx,xl
        ASM or cx,cx
        ASM je out2
        ASM mov bx,cx
        ASM shl bx,2
        ASM sub bx,4
#ifdef MR_LMM
        ASM push ds
        ASM push es
        ASM les di,DWORD PTR zg
        ASM lds si,DWORD PTR xg
#else
        ASM mov ax,ds
        ASM mov es,ax
        ASM mov di, zg
        ASM mov si, xg
#endif
        ASM add si,bx
        ASM add di,bx
        ASM mov ebx,sn
        ASM push ebp
        ASM xor ebp,ebp
    tcl2:
        ASM mov edx,ebp
        ASM lodsd
        ASM div ebx
        ASM mov ebp,edx
        ASM stosd
        ASM loop tcl2

        ASM mov eax,ebp
        ASM pop ebp
#ifdef MR_LMM
        ASM pop es
        ASM pop ds
#endif
        ASM mov sr,eax
     out2: 
        ASM cld
#endif
#if INLINE_ASM == 3
        ASM mov ecx,xl
        ASM or ecx,ecx
        ASM je out2
        ASM mov ebx,ecx
        ASM shl ebx,2
        ASM mov esi, xg
        ASM add esi,ebx
        ASM mov edi, zg
        ASM add edi,ebx
        ASM mov ebx,sn
        ASM push ebp
        ASM xor ebp,ebp
    tcl2:
        ASM sub esi,4
        ASM mov edx,ebp
        ASM mov eax,[esi]
        ASM div ebx
        ASM sub edi,4
        ASM mov ebp,edx
        ASM mov [edi],eax
        ASM dec ecx
        ASM jnz tcl2

        ASM mov eax,ebp
        ASM pop ebp