定义
若一个方程具有形式:
x
2
−
D
∗
y
2
=
1
(
D
为正整数,且
D
不是完全平方数)
x^2-D*y^2=1(D为正整数,且D不是完全平方数)
x2−D∗y2=1(D为正整数,且D不是完全平方数)
则称次方程为第一类佩尔方程,要求第一类佩尔方程的解都是正整数解(因此限制D)。
(其实有第二类,但太难了,一般用不到就不展开讲了)
性质
显然佩尔方程 x^2 - Dy^2 = 1有无数个解,化简一下就变成了
x
=
1
+
D
∗
y
2
x=\sqrt{1+D*y^2}
x=1+D∗y2
从y=1一个个试过去,如果√1+D*y^2是整数,则此时(x,y)为方程的一个解,当(x,y)最小时,也就是**最小正整数解(x1,y1),**设解为(x1,y1),(x2,y2),……。所有的解(xi,yi)可表示成如下形式:
x
i
+
√
D
∗
Y
i
=
(
x
1
+
y
1
∗
D
)
i
x^i+√D*Y^i =(x1+y1*D)^i
xi+√D∗Yi=(x1+y1∗D)i
或者由以下的递推式得到:
x
(
i
+
1
)
=
x
1
∗
x
i
+
n
∗
y
1
∗
y
i
x_(i+1)=x_1*x_i+n*y_1*y_i
x(i+1)=x1∗xi+n∗y1∗yi
y ( i + 1 ) = x 1 ∗ y i + y 1 ∗ x i y_(i+1)=x_1*y_i+y_1*x_i y(i+1)=x1∗yi+y1∗xi
由此递推式可以用前几组接就能得到后面的解
求最小正整数解的脚本如下:
import numpy as np
from collections import deque
d =
m = int(np.sqrt(d))
dq = deque()
dq.append(m)
n0 = n1 = d - m * m
m1 = m
while 1:
q, m2 = divmod(m1 + m, n1)
dq.appendleft(q)
m1 = -m2+m
n1 = (d-m1*m1)//n1
if m1 == m and n1 == n0:
break
dq.popleft()
x = 1
y = 0
for i in dq:
x1 = y + x * i
y = x
x = x1
y1=y
if x*x-d*y*y==-1:
b=(x**2+d*y**2)
y1=2*x*y
x1=b
print('x1=',x1)
print('y1=',y1)
递推式推导过程
先取方程的相邻的两个解(x1,y1),(x2,y2),显然有:
x
1
2
−
d
∗
y
1
2
=
1
x_1^2-d*y_1^2=1
x12−d∗y12=1
x
2
2
−
d
∗
y
2
2
=
1
x_2^2-d*y_2^2=1
x22−d∗y22=1
然后将两式相乘可以得到
(
x
1
2
−
d
∗
y
1
2
)
(
x
2
2
−
d
∗
y
2
2
)
=
1
(x_1^2-d*y_1^2)(x_2^2-d*y_2^2)=1
(x12−d∗y12)(x22−d∗y22)=1
再对上式进行展开
x
1
2
x
2
2
+
d
2
y
1
2
y
2
2
+
2
d
x
1
x
2
y
1
y
2
−
d
(
d
x
1
2
y
2
2
+
y
1
2
x
2
2
+
2
x
1
x
2
y
1
y
2
)
=
1
x_1^2x_2^2+d^2y_1^2y_2^2+2dx_1x2y_1y_2-d(dx_1^2y_2^2+y_1^2x_2^2+2x_1x_2y_1y_2)=1
x12x22+d2y12y22+2dx1x2y1y2−d(dx12y22+y12x22+2x1x2y1y2)=1
再合并得到
(
x
1
x
2
+
d
y
1
y
2
)
2
−
d
(
x
1
y
2
+
x
2
y
1
)
2
=
1
(x_1x_2+dy_1y_2)^2-d(x_1y_2+x_2y_1)^2=1
(x1x2+dy1y2)2−d(x1y2+x2y1)2=1
合并后发现是佩尔方程的形式,于是得到
x
3
=
x
1
x
2
+
d
y
1
y
2
x_3=x_1x_2+dy_1y_2
x3=x1x2+dy1y2
y 3 = x 1 y 2 + x 2 y 1 y_3=x_1y_2+x_2y_1 y3=x1y2+x2y1
例题
[HZNUCTF 2023 final]一步到喂
from Crypto.Util.number import *
from secret import flag,x,y
p=getPrime(512)
q=getPrime(512)
n=p*q
e=x
assert 1293023064232431070902426583269468463*pow(x,2)==105279230912868770223946474836383391725923*pow(y,2)+1293023064232431070902426583269468463
m=bytes_to_long(flag)
c=pow(m,e,n)
print(p,q)
print(c)
#9850212100620338486478343799784951721581757936621772924971218807300819701941457605399099898870264241518769370682330612103782092302148525012450902351701339
#10749429992014823019923966246896511618886613763258781706004694804949547801668777655988055847885755337127548775758133804022361510427909703124161450470578543
#66847321508502017574023222490247591875197038421108556106531660421662626233521063441647157067220450229816184622038471812597874582613385516838920632450015292570673816423432903604941781889308906893966588610214614726388822851471695742453496232748358301888465563812947038856742838097152549971517159475947566599664
仔细阅读题目条件,会发现
1293023064232431070902426583269468463*pow(x,2)==105279230912868770223946474836383391725923*pow(y,2)+1293023064232431070902426583269468463
很像佩尔方程的形式,于是我们应该会下意识的去约分,经过尝试后会发现是可行的,于是得到
x
2
−
81421
y
2
=
1
x^2-81421y^2=1
x2−81421y2=1
用上面介绍的代码可以求出整数解(x1,y1),然后求出e,最后按RSA算法的步骤求解
exp
import numpy as np
from collections import deque
import gmpy2
from Crypto.Util.number import *
d =61
m = int(np.sqrt(d))
dq = deque()
dq.append(m)
n0 = n1 = d - m * m
m1 = m
while 1:
q, m2 = divmod(m1 + m, n1)
dq.appendleft(q)
m1 = -m2+m
n1 = (d-m1*m1)//n1
if m1 == m and n1 == n0:
break
dq.popleft()
x = 1
y = 0
for i in dq:
x1 = y + x * i
y = x
x = x1
y1=y
if x*x-d*y*y==-1:
b=(x**2+d*y**2)
y1=2*x*y
x1=b
#x1=706458746417678962043621845971467865659328419354757470370053791959400836343459667183452696881749102272803495561977766862981211957547107314899837090421877848369768126470488899783907413049018564965473411482320934865623338739674557285559179500448017722749519381238449907463225595410791751764554609076679816077991393296019025846726102143252917723127850421484059625048891687413514675773883846965460650321320194376953023128886680992904874680874592398084000621265936352579085146973168432725924933271736161246147210167219273151533549827720768127049
e=x1
p=9850212100620338486478343799784951721581757936621772924971218807300819701941457605399099898870264241518769370682330612103782092302148525012450902351701339
q=10749429992014823019923966246896511618886613763258781706004694804949547801668777655988055847885755337127548775758133804022361510427909703124161450470578543
c=66847321508502017574023222490247591875197038421108556106531660421662626233521063441647157067220450229816184622038471812597874582613385516838920632450015292570673816423432903604941781889308906893966588610214614726388822851471695742453496232748358301888465563812947038856742838097152549971517159475947566599664
n=p*q
phi_n =(q-1)*(p-1)
d = gmpy2.invert(e, phi_n)
m = gmpy2.powmod(c, d, n)
flag = long_to_bytes(m)
print(flag)
#b'HZNUCTF{D0_y0u_know_p3ll_3qu4tion!!!}'