for(uint64_t i = 1; i <= log2(n); i++){
m = pow(2,i);
k_ = (p - 1)/m;
a = modExp(r,k_,p);
for(uint64_t j = 0; j < n; j+=m){
for(uint64_t k = 0; k < m/2; k++){
factor1 = result[j + k];
factor2 = modulo(modExp(a,k,p)*result[j + k + m/2],p);
result[j + k] = modulo(factor1 + factor2, p);
result[j + k+m/2] = modulo(factor1 - factor2, p);
}
}
}