Модульное возведение в степень

При обычном рекурсивном способе: a⁰=1, a ⁿ=(a ^n-1)a (mod n). Требуется n-1 операций возведения в степень по модулю т, то есть O(m²)  операций модульного умножения. Если же предварительно представить показатель п в двоичном виде  где к = [log₂ n], и, , то  где . При таком способе вычисления потребуется kмодульных возведений в квадрат и  модульных умножений. Учитывая, что каждый разряд принимает значение 0 или 1 равновероятно, можно считать, в среднем требуется 1,5 log m модульных умножений.

Алгоритм 6 Модульное возведение в степень, [Молдовян Н.А. и др., 2004].

Вход. Целое число а, 1 ≤ а < т; показатель степени п ≥ 2,

Выход. с = а^п (mod т).

1.    Положить .

2.    Для i =k-1, k-2, …, 0 выполнить следующие действия.

       2.1. Вычислить

       2.2. При n_i = 1 вычислить .

3.    Результат: с.

Код функции step (степень):

osn, st, mod – глобальные переменные основание, степень, модуль соответственно.

с = osn^st (mod mod)

void step(unsigned char c[33]){

int i,k,t,r;

int pole[256]={0};

unsigned char d[33]={0};

for(i=0;i<=31;i++){

       union{unsigned char pol;

                  struct{unsigned b0:1;

                                          unsigned b1:1;

                                          unsigned b2:1;

                                          unsigned b3:1;

                                          unsigned b4:1;

                                          unsigned b5:1;

                                          unsigned b6:1;

                                          unsigned b7:1;

                  }bit;

       }cod;

       cod.pol=st[i];

       pole[i*8+0]=cod.bit.b0;

       pole[i*8+1]=cod.bit.b1;

       pole[i*8+2]=cod.bit.b2;

       pole[i*8+3]=cod.bit.b3;

       pole[i*8+4]=cod.bit.b4;

       pole[i*8+5]=cod.bit.b5;

       pole[i*8+6]=cod.bit.b6;

       pole[i*8+7]=cod.bit.b7;

}//Представление числа в битовом формате

for(k=255;k>=0;k--){

       r=k;

       if(int(pole[k])==1)break;

}

if (int(pole[r])==0) c[0]=1;

       else

                  for(i=0;i<=31;i++){

                              c[i]=osn[i];

                  }

for(k=r-1;k>=0;k--){

       step2(c,d);

       for(t=31;t>=0;t--) c[t]=d[t];

       if(int(pole[k])==1){

       mult_mod(c,osn,d);

       for(t=31;t>=0;t--) c[t]=d[t];

       }

}

}