Кубический корень по модулю P, как мне это сделать?

Примечание добавлено позже: в алгоритме Тонелли-Шенкса и здесь предполагается, что p является простым числом. Если бы мы могли быстро вычислять модульные квадратные корни из составных модулей, мы могли бы быстро разлагать числа на множители. Прошу прощения за то, что предположил, что вы знали, что p простое число.

См. здесь или здесь. Обратите внимание, что числа по модулю p - это конечное поле с p элементами.

Изменить: см. также это (это дедушка этих документов. )

Простая часть заключается в том, что когда p = 2 mod 3, тогда все является кубом, а кубический корень из a равен просто a**((2*p-1)/3) %p.

Добавлено: Вот код, чтобы сделать все, кроме простых чисел 1 мод 9. Я постараюсь добраться до него на этих выходных. Если никто другой не доберется до него первым

#assumes p prime returns cube root of a mod p
def cuberoot(a, p):
    if p == 2:
        return a
    if p == 3:
        return a
    if (p%3) == 2:
        return pow(a,(2*p - 1)/3, p)
    if (p%9) == 4:
        root = pow(a,(2*p + 1)/9, p)
        if pow(root,3,p) == a%p:
            return root
        else:
            return None
    if (p%9) == 7:
        root = pow(a,(p + 2)/9, p)
        if pow(root,3,p) == a%p:
            return root
        else:
            return None
    else:
        print "Not implemented yet. See the second paper"

Вот полный код на чистом питоне. Если сначала рассмотреть особые случаи, то он почти так же быстр, как алгоритм Peralta.

#assumes p prime, it returns all cube roots of a mod p
def cuberoots(a, p):

    #Non-trivial solutions of x**r=1
    def onemod(p,r):
        sols=set()
        t=p-2
        while len(sols)<r:        
            g=pow(t,(p-1)//r,p)
            while g==1: t-=1; g=pow(t,(p-1)//r,p)
            sols.update({g%p,pow(g,2,p),pow(g,3,p)})
            t-=1
        return sols

    def solutions(p,r,root,a): 
        todo=onemod(p,r)
        return sorted({(h*root)%p for h in todo if pow(h*root,3,p)==a})

#---MAIN---
a=a%p

if p in [2,3] or a==0: return [a]
if p%3 == 2: return [pow(a,(2*p - 1)//3, p)] #One solution

#There are three or no solutions 

#No solution
if pow(a,(p-1)//3,p)>1: return []

if p%9 == 7:                                #[7, 43, 61, 79, 97, 151]   
    root = pow(a,(p + 2)//9, p)
    if pow(root,3,p) == a: return solutions(p,3,root,a)
    else: return []

if p%9 == 4:                                #[13, 31, 67, 103, 139]
    root = pow(a,(2*p + 1)//9, p) 
    print(root)
    if pow(root,3,p) == a: return solutions(p,3,root,a)        
    else: return []        
            
if p%27 == 19:                              #[19, 73, 127, 181]
    root = pow(a,(p + 8)//27, p)
    return solutions(p,9,root,a)

if p%27 == 10:                              #[37, 199, 307]
    root = pow(a,(2*p +7)//27, p)  
    return solutions(p,9,root,a) 

#We need a solution for the remaining cases
return tonelli3(a,p,True)

Расширение алгоритма Тонелли-Шэнка.

def tonelli3(a,p,many=False):

    def solution(p,root):
        g=p-2
        while pow(g,(p-1)//3,p)==1: g-=1  #Non-trivial solution of x**3=1
        g=pow(g,(p-1)//3,p)
        return sorted([root%p,(root*g)%p,(root*g**2)%p])

#---MAIN---
a=a%p
if p in [2,3] or a==0: return [a]
if p%3 == 2: return [pow(a,(2*p - 1)//3, p)] #One solution

#No solution
if pow(a,(p-1)//3,p)>1: return []

#p-1=3**s*t
s=0
t=p-1
while t%3==0: s+=1; t//=3
 
#Cubic nonresidu b
b=p-2
while pow(b,(p-1)//3,p)==1: b-=1

c,r=pow(b,t,p),pow(a,t,p)    
c1,h=pow(c,3**(s-1),p),1    
c=pow(c,p-2,p) #c=inverse modulo p

for i in range(1,s):
    d=pow(r,3**(s-i-1),p)
    if d==c1: h,r=h*c,r*pow(c,3,p)
    elif d!=1: h,r=h*pow(c,2,p),r*pow(c,6,p)           
    c=pow(c,3,p)
    
if (t-1)%3==0: k=(t-1)//3
else: k=(t+1)//3

r=pow(a,k,p)*h
if (t-1)%3==0: r=pow(r,p-2,p) #r=inverse modulo p

if pow(r,3,p)==a: 
    if many: 
        return solution(p,r)
    else: return [r]
else: return [] 

Вы можете протестировать его, используя:

test=[(17,1459),(17,1000003),(17,10000019),(17,1839598566765178548164758165715596714561757494507845814465617175875455789047)]

for a,p in test:
    print "y^3=%s modulo %s"%(a,p)
    sol=cuberoots(a,p)
    print "p%s3=%s"%("%",p%3),sol,"--->",map(lambda t: t^3%p,sol)

который должен дать (быстро):

y^3=17 по модулю 1459
p%3=1 [483, 329, 647] ---› [17, 17, 17]
y^3=17 по модулю 1000003
p%3= 1 [785686, 765339, 448981] ---› [17, 17, 17]
y^3=17 по модулю 10000019
p%3=2 [5188997] ---› [17]
у ^ 3 = 17 по модулю 1839598566765178548164758165715596714561757494507845814465617175875455789047
р% 3 = 1 [753801617033579226225229608063663938352746555486783903392457865386777137044, 655108821219252496141403783945148550782812009720868259303598196387356108990, 430688128512346825798124773706784225426198929300193651769561114101322543013] ---> [17, 17, 17]