3-2.針對大的\(r\)，分解RSA公鑰\(n=p^rq\)
方法 | 範例 |
問題敘述 |
定理
令\(N=p^rq\)，對某些的\(c\)，\(q<p^c\)。
若已知\(N,r,c\)由執行時間為\(\displaystyle exp\left(\frac{c+1}{r+c}\cdot logp\right)\cdot O(\gamma)\)的演算法能分解出因數\(p\)。
其中\(\gamma\)是矩陣維度\(O(r^2)\)，矩陣元素大小為\(O(\gamma logN)\)執行LLL所需的時間。
這個演算法是決定性的和需要多項式的空間。
[證明]
\(N=p^rq\)，對某些的\(c\)，\(q<p^c\)
令\(d=2r(r+c)\)，由引理1可知整數\(P\)滿足\(|\; P-p|\;<p^{1-\frac{c+1}{r+c}}\)
就能在\(O((logN)^2d^4)\)的時間內找到\(N\)的分解。 | \(197213=61^2 \cdot 53\)要找出\(197213\)的因數\(61\)
\(N=197213\)，\(r=2\)，\(p=61\)，\(q=53\)
假設\(c=1\)，\(53<61^c\)
矩陣維度\(d=2r(r+c)=2\cdot 2(2+1)=12\)
論文範例取較小的\(d\)值，\(d=9\) |
步驟1：列舉因數\(p\)的近似值\(P\) |
令\(\displaystyle ε=\frac{c+1}{r+c}\)
For \(k=1,2,\ldots,\frac{log_2N}{r}\) do
For \(j=0,1,\ldots,2^{ε}\) do
(設\(P=2^k+j\cdot 2^{(1-ε)k}\)
使用近似的\(P\)，執行引理1的演算法
) | \( \displaystyle ε=\frac{c+1}{r+c}=\frac{1+1}{2+1}=\frac{2}{3}\)
\(\displaystyle k=5,j=9,P=2^5+9*2^{(1-2/3)5}=59\)
取\(P=59\)當作因數\(p\)的近似值 |
引理1：
給定\(N=p^r q\)和假設\(q<p^c\)對某些\(c\)，更進一步假設\(P\)是符合\(\displaystyle |\; P-p|\;<p^{1-\frac{c}{r+c}-2\frac{r}{d}}\)的整數
從\(N,r,c\)和\(P\)使用維度\(d\)的lattice執行LLL方法，能計算出因數\(p\)。 |
步驟2.設同餘方程式，計算參數\(X\) |
設同餘方程式為\(f(x)=(P+x)^r\equiv 0 \pmod{p^r}\)
且\(p(x)\)為monic(最高次方項係數為1)。
利用LLL方法可以找出比邊界\(X\)還小的解\(x_0=p-P\)(\(\displaystyle |\;x_0|\;<X\))
使得\(f(x_0)\equiv 0 \pmod{p^r}\)，進而得到公鑰\(N\)的一個因數\(p=P+x_0\)。 | \(f(x)=(35+x)^2 \equiv 0\pmod{p^2}\)
\(\displaystyle X=p^{1-\frac{c+1}{r+c}}=p^{1-ε}\)
因為\(p\)未知，假設\(p,q\)相近，\(\displaystyle N=p^rq\approx p^{r+1}\)，\(\displaystyle p\approx N^{\frac{1}{r+1}}\)
\(\displaystyle X\approx N^{\frac{1-ε}{r+1}}= N^{\frac{1-2/3}{2+1}}=197213^{\frac{1}{9}}\approx 3.875\)
取\(X=3\) |
步驟3.產生矩陣\(M\) |
\(g_{i,k}(xX)=N^{m-k}(xX)^i f(xX)^k\)
\(g_{i,k}(xX)\)，\(i=0,1,\ldots,r-1\)和\(k=0,1,\ldots,m-1\)
\(g_{j,m}(xX)\)，\(j=0,1,\ldots,d-mr-1\)
注意到對所有\(i,k\)，\(g_{i,k}(x_0)\equiv 0\pmod{p^{rm}}\)
[證明]
\(g_{i,k}(x_0)=N^{m-k}(x_0X)^i(P+x_0)^{rk}\)
\(=(p^rq)^{m-k}(x_0X)^i(p)^{rk}\)
\(=p^{rm}(q^{m-k}(x_0X)^i)\)
\(\equiv 0\pmod{p^{rm}}\)
注意到對所有\(j,m\)，\(g_{j,m}(x_0)\equiv 0\pmod{p^{rm}}\)
[證明]
\(g_{j,m}(x_0)=N^{m-m}(x_0X)^j(P+x_0)^{rm}\)
\(=(x_0X)^jp^{rm}\)
\(\equiv 0\pmod{p^{rm}}\)
將\(g_{i,k}(xX),g_{j,m}(xX)\)多項式依\(x\)次方形成係數矩陣\(M\)。
令\(M\)是\(L\)的基底向量所形成的矩陣，注意到\(M\)為下三角矩陣，所以\(L\)的行列式值為\(M\)的對角線元素相乘，得到
\(\displaystyle det(M)=\left(\prod_{k=0}^{m-1}\prod_{i=0}^{r-1}N^{m-k}\right)\left(\prod_{j=0}^{d-1}X^j\right)\)
\(\displaystyle =\left(\prod_{i=0}^{r-1}N^m\cdot N^{m-1}\ldots N^1\right) \cdot (X^0 \cdot X^1 \ldots X^{d-1})\)
\(\displaystyle =\left(\prod_{i=0}^{r-1}N^{m(m+1)/2}\right)X^{d(d-1)/2}\)
\(\displaystyle =N^{rm(m+1)/2}X^{d(d-1)/2}\)
\(\displaystyle <N^{rm(m+1)/2}X^{d^2/2}\) | \(r=2\)
\(\displaystyle m=\lfloor\; \frac{d}{r+c}-\frac{1}{2} \rfloor\;=\lfloor\; \frac{12}{2+1}-\frac{1}{2} \rfloor\;=3\)
\(d=9\)
產生\(g_{i,k}(xX)\)多項式
\(i=0,k=0,g_{0,0}(xX)=N^3(xX)^0f(xX)^0=\)\(\bbox[border:1px solid black]{N^3}\)
\(i=1,k=0,g_{1,0}(xX)=N^3(xX)^1f(xX)^0=\)\(\bbox[border:1px solid black]{N^3X}x\)
\(i=0,k=1,g_{0,1}(xX)=N^2(xX)^0f(xX)^1=N^2(P+Xx)^2\)
\(=N^2P^2+2N^2PXx+\)\(\bbox[border:1px solid black]{N^2X^2}x^2\)
\(i=1,k=1,g_{1,1}(xX)=N^2(xX)^1f(xX)^1=N^2Xx(P+Xx)^2\)
\(=N^2P^2Xx+2N^2PX^2x^2+\bbox[border:1px solid black]{N^2X^3}x^3\)
\(i=0,k=2,g_{0,2}(xX)=N(xX)^0f(xX)^2=N(P+Xx)^4\)
\(=NP^4+4NP^3Xx+6NP^2X^2x^2+4NPX^3x^3+\bbox[border:1px solid black]{NX^4}x^4\)
\(i=1,k=2,g_{1,2}(xX)=N(xX)^1f(xX)^2=NXx(P+Xx)^4\)
\(=NP^4Xx+4NP^3X^2x^2+6NP^2X^3x^3+4NPX^4x^4+\bbox[border:1px solid black]{NX^5}x^5\)
產生\(g_{j,m}(xX)\)多項式
\(j=0,g_{0,3}(xX)=(xX)^0f(xX)^3=(P+Xx)^6\)
\(=P^6+6P^5Xx+15P^4X^2x^2+20P^3X^3x^3+15P^2X^4x^4+6PX^5x^5+\bbox[border:1px solid black]{X^6}x^6\)
\(j=1,g_{1,3}(xX)=(xX)f(xX)^3=Xx(P+Xx)^6\)
\(=P^6Xx+6P^5X^2x^2+15P^4X^3x^3+20P^3X^4x^4+15P^2X^5x^5+6PX^6x^6+\bbox[border:1px solid black]{X^7}x^7\)
\(j=2,g_{2,3}(xX)=(xX)^2f(xX)^3=X^2x^2(P+Xx)^6\)
\(=P^6X^2x^2+6P^5X^3x^3+15P^4X^4x^4+20P^3X^5x^5+15P^2X^6x^6+6PX^7x^7+\bbox[border:1px solid black]{X^8}x^8\) |
\(M=\matrix{&\matrix{1 &x &x^2 &x^3& x^4 & x^5& x^6& x^7& x^8}\cr
\matrix{g_{0,0}(xX)\cr g_{1,0}(xX)\cr g_{0,1}(xX)\cr g_{1,1}(xX)\cr g_{0,2}(xX)\cr g_{1,2}(xX)\cr g_{0,3}(xX)\cr g_{1,3}(xX)\cr g_{2,3}(xX)}&\left[\matrix{N^3&&&&&&&&\cr
0&XN^3&&&&&&&\cr
*&*&X^2N^2&&&&&&\cr
0&*&*&X^3N^2&&&&&\cr
*&*&*&*&X^4N&&&&\cr
0&*&*&*&*&X^5N&&&\cr
*&*&*&*&*&*&X^6&&\cr
0&*&*&*&*&*&*&X^7&\cr
0&0&*&*&*&*&*&*&X^8}\right]}\)
*代表非零數字
\(det(M)=N^3 \cdot XN^3 \cdot X^2N^2\cdot X^3N^2\cdot X^4N\cdot X^5N\cdot X^6\cdot X^7\cdot X^8\)
\(=N^{2\cdot(3+2+1)}X^{0+1+2+\ldots+7+8}\)
\(=N^{12}X^{36}\) |
步驟4.經LLL化簡後的短向量產生不需要同餘\(p^{rm}\)的方程式，得到公鑰\(N\)的因數\(p\) |
矩陣\(M\)經LLL化簡為\(B\)
\(B=LLL(M)\)
lattice經LLL化簡後第一行\(b_1\)為整個lattice中較短向量
所形成的方程式不需要再同餘\(p^{rm}\)
－－－－－－－－－－－－－－－－－－－
令\(h(x)\in Z[x]\)是次方為\(n\)的多項式，假設
(1)\(h(x_0)\equiv 0\pmod{p^{rm}}\)，\(r,m\)為正整數，\(|\;x_0|\;<X\)
(2)\(\displaystyle \Vert\;h(xX)\Vert\;<\frac{p^{rm}}{\sqrt{d}}\)
則在整數下\(h(x_0)=0\)。
[證明]
\(\displaystyle |\;h(x_0)|\;=|\;\sum a_ix_0^i|\;=|\;\sum a_i X^i \left(\frac{x_0}{X}\right)^i|\;\)
\(\displaystyle \le \sum |\;a_iX^i \left(\frac{x_0}{X}\right)^i|\;\le \sum |\;a_i X^i|\;\)
\(\le \sqrt{d} \Vert\;h(xX)\Vert\;<p^{rm}\)
因為\(h(x_0)\equiv 0\pmod{p^{rm}}\)，所以\(h(x_0)=0\)
－－－－－－－－－－－－－－－－－－－
向量\(b_1\)是由多項式\(g_{i,k}(xX)\)係數的線性組合而成\(h(xX)\)，\(\Vert\;b_1 \Vert\;=\Vert\;h(xX) \Vert\;\)
\(h(x)\)也會是多項式\(g_{i,k}(x)\)的線性組合，因為\(g_{i,k}(x_0)\equiv 0\pmod{p^{rm}}\)，所以\(h(x_0)\equiv 0\pmod{p^{rm}}\)
假設lattice \(L\)經LLL化簡後向量\(b_1,b_2,\ldots,b_n\)，則\(\Vert\;b_1 \Vert\;\le 2^{(n-1)/4}\cdot(det(L))^{1/n}\)。
證明在https://math.pro/db/viewthread.php?tid=3498&page=1#pid23067
論文將\(2^{(n-1)/4}\)改成\(2^{d/2}\)，\(d\)為lattice矩陣維度
\(\Vert\;b_1 \Vert\;\le 2^{d/2}\cdot(det(L))^{1/d}\)
\(\Vert\;b_1 \Vert\;^d\le 2^{d^2/2}\cdot det(L)\)
將行列式值\( det(L)<N^{rm(m+1)/2}X^{d^2/2}\)代入得到
\(\Vert\;b_1 \Vert\;^d\le 2^{d^2/2}N^{rm(m+1)/2}X^{d^2/2}\)
向量\(b_1\)是某個\(h(xX)\)多項式的係數向量，滿足\(\Vert\;b_1 \Vert\;=\Vert\;h(xX) \Vert\;\)
更進一步，因為\(h(xX)\)是\(g_{i,k}(xX)\)多項式的線性組合，
所以\(h(x)\)也是\(g_{i,k}(x)\)多項式的線性組合
將\(\Vert\;b_1 \Vert\;=\Vert\;h(xX) \Vert\;\)代入
\(\Vert\;h(xX) \Vert\;^d\le 2^{d^2/2}\cdot N^{rm(m+1)/2}X^{d^2/2}\)
若\(\displaystyle \Vert\;h(xX)\Vert\;<\frac{p^{rm}}{\sqrt{d}}\)，分母的\(\sqrt{d}\)在接下來計算中影響較小，為了簡化就忽略了。
\(\displaystyle \Vert\;h(xX)\Vert\;^d\le 2^{d^2/2}\cdot N^{rm(m+1)/2}X^{d^2/2}<p^{rmd}\)
\((2X)^{d^2/2}<p^{rmd}N^{-rm(m+1)/2}\)
假設對某個\(c\)，\(q<p^c\)，\(N=p^rq<p^{r+c}\)
\((2X)^{d^2/2}<p^{rmd}p^{-r(r+c)m(m+1)/2}\)
\((2X)^{d^2/2}<p^{rmd-r(r+c)m(m+1)/2}\)
為了讓\(X\)值越大越好，可以取較弱的近似\(P\)值
利用配方法求\(p\)的最大次方\(rmd-r(r+c)m(m+1)/2\)
\(\displaystyle =rmd-\frac{r}{2}(r+c)m(m+1)\)
\(\displaystyle =-\frac{r}{2}\left[(r+c)m^2+(r+c-2d)m\right]\)
\(\displaystyle =-\frac{r}{2}(r+c)\left(m^2+\frac{r+c-2d}{r+c}m\right)\)
\(\displaystyle =-\frac{r}{2}(r+c)\left(m+\frac{r+c-2d}{2(r+c)}\right)^2+\frac{r}{2}(r+c)\left(\frac{r+c-2d}{2(r+c)}\right)^2\)
\(\displaystyle =\frac{r}{2}(r+c)\left(m+\frac{1}{2}-\frac{d}{r+c}\right)^2+\frac{r}{2}(r+c)\left(\frac{r+c-2d}{2(r+c)}\right)^2\)
\(m\)的最佳值在\(\displaystyle m_0=\lfloor\; \frac{d}{r+c}-\frac{1}{2} \rfloor\;\)
我們也許能選\(d_0\)使得\(\displaystyle \frac{d_0}{r+c}-\frac{1}{2}\)在\(\displaystyle \frac{1}{2r+c}\)的整數範圍內。
將\(m=m_0\)和\(d=d_0\)代入和繁瑣的計算後得到
\(\displaystyle X<\frac{1}{2}p^{1-\frac{c}{r+c}-\frac{r}{d}(1+\delta)}\)，其中\(\displaystyle \delta=\frac{1}{r+c}-\frac{r+c}{4d}\)
因為\(\delta<1\)我們得到稍微弱但更吸引人的界限
\(\displaystyle X<p^{1-\frac{c}{r+c}-2\frac{r}{d}}\)
令\(d=2r(r+c)\)
\(\displaystyle X<p^{1-\frac{c}{r+c}-2\frac{r}{2r(r+c)}}\)
\(\displaystyle X<p^{1-\frac{c+1}{r+c}}\) | lattice經LLL化簡後第一行\(b_1\)為整個lattice中較短向量
\(b_1=[83412628,-91271724,-27587799,-100623168,\)
\([-100709649,89252928,78155361,202680225,228782070]\)
產生不需要同餘\(p^{rm}\)的方程式
\(\displaystyle h(x)=83412628-91271724 \left(\frac{x}{3}\right)-27587799 \left(\frac{x}{3}\right)^2\)
\(\displaystyle-100623168 \left(\frac{x}{3}\right)^3-100709649 \left(\frac{x}{3}\right)^4+89252928 \left(\frac{x}{3}\right)^5\)
\(\displaystyle+78155361\left(\frac{x}{3}\right)^6+202680225\left(\frac{x}{3}\right)^7+228782070 \left(\frac{x}{3}\right)^8\)
\(=83412628-30423908x-3065311x^2-3726784x^3\)
\(-1243329x^4+367296x^5+107209x^6+92675x^7+34870x^8\)
\(=(-2+x)^2(20853157+13247180x+7267563x^2+\)
\(3024072x^3+896349x^4+232155x^5+34870x^6)\)
得到整數解\(x=2\)
得到因數\(p=P+x=59+2=61\)
\(N\)可分解成\(197213=61^2 \cdot 53\) |
參考資料：
Boneh, D., Durfee, G., Howgrave-Graham, N. (1999). Factoring N = p r q for Large r . In: Wiener, M. (eds) Advances in Cryptology — CRYPTO’ 99. CRYPTO 1999. Lecture Notes in Computer Science, vol 1666. Springer, Berlin, Heidelberg.
https://doi.org/10.1007/3-540-48405-1_21
https://link.springer.com/chapte ... 5-1_21#chapter-info
請下載LLL.zip，解壓縮後將LLL.mac放到C:\maxima-5.46.0\share\maxima\5.46.0\share目錄下
要先載入LLL.mac才能使用LLL指令
load("LLL.mac");
C:/maxima-5.47.0/share/maxima/5.47.0/share/LLL.mac
要因數分解的公鑰\(N\)
N:197213;
\(197213\)
\(N=p^r q\)，因數\(p\)的次方\(r\)
r:2;
\(2\)
假設\(q<p^c\)
c:1;
\(1\)
矩陣維度
d:2*r*(r+c);
\(12\)
參數\(m\)
m:floor(d/(r+c)-1/2);
\(3\)
為了論文範例將矩陣維度改為9
d:9;
\(9\)
列舉\(P\)的ε
epsilon: (c+1)/(r+c);
\(\displaystyle \frac{2}{3}\)
列舉因數\(p\)的近似值\(P\)
for k:1 thru floor(log(N)/(log(2)*r)) do
(steps:floor(float(2^(epsilon*k))),
for j:0 thru steps do
(P:2^k+j*floor(float(2^((1-epsilon)*k))),
print("k=",k,",j=",j,",P=2"^k,"+",j,"*2"^((1-epsilon)*k),"=",P)
)
);
\(k=1,j=0,P=2+0*2^{1/3}=2\)
\(k=1,j=1,P=2+1*2^{1/3}=3\)
\(k=2,j=0,P=2^2+0*2^{2/3}=4\)
\(k=2,j=1,P=2^2+1*2^{2/3}=5\)
\(k=2,j=2,P=2^2+2*2^{2/3}=6\)
\(k=3,j=0,P=2^3+0*2=8\)
\(k=3,j=1,P=2^3+1*2=10\)
\(k=3,j=2,P=2^3+2*2=12\)
\(k=3,j=3,P=2^3+3*2=14\)
\(k=3,j=4,P=2^3+4*2=16\)
\(k=4,j=0,P=2^4+0*2^{4/3}=16\)
\(k=4,j=1,P=2^4+1*2^{4/3}=18\)
\(k=4,j=2,P=2^4+2*2^{4/3}=20\)
\(k=4,j=3,P=2^4+3*2^{4/3}=22\)
\(k=4,j=4,P=2^4+4*2^{4/3}=24\)
\(k=4,j=5,P=2^4+5*2^{4/3}=26\)
\(k=4,j=6,P=2^4+6*2^{4/3}=28\)
\(k=5,j=0,P=2^5+0*2^{5/3}=32\)
\(k=5,j=1,P=2^5+1*2^{5/3}=35\)
\(k=5,j=2,P=2^5+2*2^{5/3}=38\)
\(k=5,j=3,P=2^5+3*2^{5/3}=41\)
\(k=5,j=4,P=2^5+4*2^{5/3}=44\)
\(k=5,j=5,P=2^5+5*2^{5/3}=47\)
\(k=5,j=6,P=2^5+6*2^{5/3}=50\)
\(k=5,j=7,P=2^5+7*2^{5/3}=53\)
\(k=5,j=8,P=2^5+8*2^{5/3}=56\)
\(\bbox[border:1px solid black]{k=5,j=9,P=2^5+9*2^{5/3}=59}\)
\(k=5,j=10,P=2^5+10*2^{5/3}=62\)
\(k=6,j=0,P=2^6+0*2^2=64\)
\(k=6,j=1,P=2^6+1*2^2=68\)
\(k=6,j=2,P=2^6+2*2^2=72\)
\(k=6,j=3,P=2^6+3*2^2=76\)
\(k=6,j=4,P=2^6+4*2^2=80\)
\(k=6,j=5,P=2^6+5*2^2=84\)
\(k=6,j=6,P=2^6+6*2^2=88\)
\(k=6,j=7,P=2^6+7*2^2=92\)
\(k=6,j=8,P=2^6+8*2^2=96\)
\(k=6,j=9,P=2^6+9*2^2=100\)
\(k=6,j=10,P=2^6+10*2^2=104\)
\(k=6,j=11,P=2^6+11*2^2=108\)
\(k=6,j=12,P=2^6+12*2^2=112\)
\(k=6,j=13,P=2^6+13*2^2=116\)
\(k=6,j=14,P=2^6+14*2^2=120\)
\(k=6,j=15,P=2^6+15*2^2=124\)
\(k=6,j=16,P=2^6+16*2^2=128\)
\(k=7,j=0,P=2^7+0*2^{7/3}=128\)
\(k=7,j=1,P=2^7+1*2^{7/3}=133\)
\(k=7,j=2,P=2^7+2*2^{7/3}=138\)
\(k=7,j=3,P=2^7+3*2^{7/3}=143\)
\(k=7,j=4,P=2^7+4*2^{7/3}=148\)
\(k=7,j=5,P=2^7+5*2^{7/3}=153\)
\(k=7,j=6,P=2^7+6*2^{7/3}=158\)
\(k=7,j=7,P=2^7+7*2^{7/3}=163\)
\(k=7,j=8,P=2^7+8*2^{7/3}=168\)
\(k=7,j=9,P=2^7+9*2^{7/3}=173\)
\(k=7,j=10,P=2^7+10*2^{7/3}=178\)
\(k=7,j=11,P=2^7+11*2^{7/3}=183\)
\(k=7,j=12,P=2^7+12*2^{7/3}=188\)
\(k=7,j=13,P=2^7+13*2^{7/3}=193\)
\(k=7,j=14,P=2^7+14*2^{7/3}=198\)
\(k=7,j=15,P=2^7+15*2^{7/3}=203\)
\(k=7,j=16,P=2^7+16*2^{7/3}=208\)
\(k=7,j=17,P=2^7+17*2^{7/3}=213\)
\(k=7,j=18,P=2^7+18*2^{7/3}=218\)
\(k=7,j=19,P=2^7+19*2^{7/3}=223\)
\(k=7,j=20,P=2^7+20*2^{7/3}=228\)
\(k=7,j=21,P=2^7+21*2^{7/3}=233\)
\(k=7,j=22,P=2^7+22*2^{7/3}=238\)
\(k=7,j=23,P=2^7+23*2^{7/3}=243\)
\(k=7,j=24,P=2^7+24*2^{7/3}=248\)
\(k=7,j=25,P=2^7+25*2^{7/3}=253\)
\(k=8,j=0,P=2^8+0*2^{8/3}=256\)
\(k=8,j=1,P=2^8+1*2^{8/3}=262\)
\(k=8,j=2,P=2^8+2*2^{8/3}=268\)
\(k=8,j=3,P=2^8+3*2^{8/3}=274\)
\(k=8,j=4,P=2^8+4*2^{8/3}=280\)
\(k=8,j=5,P=2^8+5*2^{8/3}=286\)
\(k=8,j=6,P=2^8+6*2^{8/3}=292\)
\(k=8,j=7,P=2^8+7*2^{8/3}=298\)
\(k=8,j=8,P=2^8+8*2^{8/3}=304\)
\(k=8,j=9,P=2^8+9*2^{8/3}=310\)
\(k=8,j=10,P=2^8+10*2^{8/3}=316\)
\(k=8,j=11,P=2^8+11*2^{8/3}=322\)
\(k=8,j=12,P=2^8+12*2^{8/3}=328\)
\(k=8,j=13,P=2^8+13*2^{8/3}=334\)
\(k=8,j=14,P=2^8+14*2^{8/3}=340\)
\(k=8,j=15,P=2^8+15*2^{8/3}=346\)
\(k=8,j=16,P=2^8+16*2^{8/3}=352\)
\(k=8,j=17,P=2^8+17*2^{8/3}=358\)
\(k=8,j=18,P=2^8+18*2^{8/3}=364\)
\(k=8,j=19,P=2^8+19*2^{8/3}=370\)
\(k=8,j=20,P=2^8+20*2^{8/3}=376\)
\(k=8,j=21,P=2^8+21*2^{8/3}=382\)
\(k=8,j=22,P=2^8+22*2^{8/3}=388\)
\(k=8,j=23,P=2^8+23*2^{8/3}=394\)
\(k=8,j=24,P=2^8+24*2^{8/3}=400\)
\(k=8,j=25,P=2^8+25*2^{8/3}=406\)
\(k=8,j=26,P=2^8+26*2^{8/3}=412\)
\(k=8,j=27,P=2^8+27*2^{8/3}=418\)
\(k=8,j=28,P=2^8+28*2^{8/3}=424\)
\(k=8,j=29,P=2^8+29*2^{8/3}=430\)
\(k=8,j=30,P=2^8+30*2^{8/3}=436\)
\(k=8,j=31,P=2^8+31*2^{8/3}=442\)
\(k=8,j=32,P=2^8+32*2^{8/3}=448\)
\(k=8,j=33,P=2^8+33*2^{8/3}=454\)
\(k=8,j=34,P=2^8+34*2^{8/3}=460\)
\(k=8,j=35,P=2^8+35*2^{8/3}=466\)
\(k=8,j=36,P=2^8+36*2^{8/3}=472\)
\(k=8,j=37,P=2^8+37*2^{8/3}=478\)
\(k=8,j=38,P=2^8+38*2^{8/3}=484\)
\(k=8,j=39,P=2^8+39*2^{8/3}=490\)
\(k=8,j=40,P=2^8+40*2^{8/3}=496\)
done
要試驗很多次才能找到合適的近似值\(P\)，為了簡化過程選擇\(P=59\)
P:59;
\(59\)
希望能找到\(|\;x|\;<X\)，\(f(x)\equiv 0\pmod{p^r}\)
X:floor(N^((1-epsilon)/(r+1)));
\(3\)
同餘方程式\(f(x)=(P+x)^r\equiv 0\pmod{p^r}\)
fx: ('P+x)^r;
\((x+P)^2\)
將\(x\)以\(Xx\)代替，\(f(Xx)=(P+Xx)^r \equiv 0\pmod{p^r}\)
fXx:subst(x=x*'X,fx);
\((Xx+P)^2\)
以\(x\)升冪排序顯示
powerdisp:true;
\(true\)
\(g(xX)\)多項式
gxX:[];
[]
產生\(g_{i,k}(xX)=N^{m-k}(Xx)^if^k(xX)\)多項式，\(i=0,\ldots,r-1\)，\(k=0,\ldots,m-1\)
for k:0 thru m-1 do
(for i:0 thru r-1 do
(print("i=",i,",k=",k,",g",i,",",k,"(xX)=N"^(m-k),"*","(xX)"^i,"*","f(xX)"^k,"=",gik:'N^(m-k)*(x*'X)^i*fXx^k,"=",expand(gik)),
gxX:append(gxX,[gik])
)
);
\(i=0,k=0,g_{0,0}(xX)=N^3*1*1=N^3=N^3\)
\(i=1,k=0,g_{1,0}(xX)=N^3*(xX)*1=N^3Xx=N^3Xx\)
\(i=0,k=1,g_{0,1}(xX)=N^2*1*f(xX)=N^2(P+Xx)^2=N^2P^2+2N^2PXx+N^2X^2x^2\)
\(i=1,k=1,g_{1,1}(xX)=N^2*(xX)*f(xX)=N^2Xx(P+Xx)^2=N^2P^2Xx+2N^2PX^2x^2+N^2X^3x^3\)
\(i=0,k=2,g_{0,2}(xX)=N*1*f(xX)^2=N(P+Xx)^4=NP^4+4NP^3Xx+6NP^2X^2x^2+4NPX^3x^3+NX^4x^4\)
\(i=1,k=2,g_{1,2}(xX)=N*(xX)*f(xX)^2=NXx(P+Xx)^4=NP^4Xx+4NP^3X^2x^2+6NP^2X^3x^3+4NPX^4x^4+NX^5x^5\)
done
產生\(g_{j,m}=x^jf^m(xX)\)多項式，\(j=0,\ldots,d-mr-1\)
for j:0 thru d-m*r-1 do
(print("j=",j,",g",j,",",m,"(xX)=","(xX)"^j,"*f(xX)"^m,"=",gim: (x*'X)^j*fXx^m,"=",expand(gim)),
gxX:append(gxX,[gim])
);
\(j=0,g_{0,3}(xX)=1*f(xX)^3=(P+Xx)^6=P^6+6P^5Xx+15P^4X^2x^2+20P^3X^3x^3+15P^2X^4x^4+6PX^5x^5+X^6x^6\)
\(j=1,g_{1,3}(xX)=(xX)*f(xX)^3=Xx(P+Xx)^6=P^6Xx+6P^5X^2x^2+15P^4X^3x^3+20P^3X^4x^4+15P^2X^5x^5+6PX^6x^6+X^7x^7\)
\(j=2,g_{2,3}(xX)=(xX)^2*f(xX)^3=X^2x^2(P+Xx)^6=P^6X^2x^2+6P^5X^3x^3+15P^4X^4x^4+20P^3X^5x^5+15P^2X^6x^6+6PX^7x^7+X^8x^8\)
done
全部的\(g(xX)\)多項式
gxX;
\([N^3,N^3Xx,N^2(P+Xx)^2,N^2Xx(P+Xx)^2,N(P+Xx)^4,NXx(P+Xx)^4,(P+Xx)^6,Xx(P+Xx)^6,X^2x^2(P+Xx)^6]\)
\(x^1,\ldots,x^{d-1}\)
xpower:create_list(x^i,i,1,d-1);
\([x,x^2,x^3,x^4,x^5,x^6,x^7,x^8]\)
取\(g(xX)\)多項式係數(常數項在最後一行)
M:augcoefmatrix(gxX,xpower);
\(\left[\matrix{0&0&0&0&0&0&0&0&N^3\cr
N^3X&0&0&0&0&0&0&0&0\cr
2N^2PX&N^2X^2&0&0&0&0&0&0&N^2P^2\cr
N^2P^2X&2N^2PX^2&N^2X^3&0&0&0&0&0&0\cr
4NP^3X&6NP^2X^2&4NPX^3&NX^4&0&0&0&0&NP^4\cr
NP^4X&4NP^3X^2&6NP^2X^3&4NPX^4&NX^5&0&0&0&0\cr
6P^5X&15P^4X^2&20P^3X^3&15P^2X^4&6PX^5&X^6&0&0& P^6\cr
P^6X&6P^5X^2&15P^4X^3&20P^3X^4&15P^2X^5&6P X^6&X^7&0&0\cr
0& P^6X^2&6P^5X^3&15P^4X^4&20P^3X^5&15P^2X^6&6PX^7&X^8&0}\right]\)
將常數項移到第一行
M:addcol(col(M,d),submatrix(M,d));
\(\left[\matrix{N^3&0&0&0&0&0&0&0&0\cr
0&N^3X&0&0&0&0&0&0&0\cr
N^2P^2&2N^2PX&N^2X^2&0&0&0&0&0&0\cr
0&N^2P^2X&2N^2PX^2&N^2X^3&0&0&0&0&0\cr
NP^4&4NP^3X&6NP^2X^2&4NPX^3&NX^4&0&0&0&0\cr
0&NP^4X&4NP^3X^2&6NP^2X^3&4NPX^4&NX^5&0&0&0\cr
P^6&6P^5X&15P^4X^2&20P^3X^3&15P^2X^4&6PX^5&X^6&0&0\cr
0& P^6X&6P^5X^2&15P^4X^3&20P^3X^4&15P^2X^5&6PX^6&X^7&0\cr
0&0& P^6X^2&6P^5X^3&15P^4X^4&20P^3X^5&15P^2X^6&6PX^7&X^8}\right]\)
將\(N=197213,X=3\)代入矩陣\(M\)
M:ev(M,[N=N,P=P,X=X]);
\(\left[\matrix{7670198773742597&0&0&0&0&0&0&0&0\cr
0&23010596321227791&0&0&0&0&0&0&0\cr
135386419411489&13768110448626&350036706321&0&0&0&0&0&0\cr
0&406159258234467&41304331345878&1050110118963&0&0&0&0&0\cr
2389701114893&486040904724&37070916462&1256641236&15974253&0&0&0&0\cr
0&7169103344679&1458122714172&111212749386&3769923708&47922759&0&0&0\cr
42180533641&12868637382&1635843735&110904660&4229415&86022&729&0&0\cr
0&126541600923&38605912146&4907531205&332713980&12688245&258066&2187&0\cr
0&0&379624802769&115817736438&14722593615&998141940&38064735&774198&6561}\right]\)
LLL化簡
B: LLL(M);
\(\left[\matrix{83412628&-91271724&-27587799&-100623168&-100709649&89252928&78155361&202680225&228782070\cr
164609528&-217976916&-267622146&110184705&219274371&225169146&-62478216&-7989111&-216513\cr
-78834856&97014168&-7341966&273951288&-166154733&-57299643&-299831139&-30143421&120564936\cr
-148796076&254769516&149231061&-295864218&87910272&28543023&-249134292&82347111&-96958458\cr
1763788&20402880&-16597377&-47384811&-126043452&287403633&82562166&-398327058&234910044\cr
23658984&-157701924&142819218&-8117631&410281605&-146082852&-373406193&-393869952&369922302\cr
-92491784&338523084&-376090650&50862033&300067821&-401077818&111970026&-180821160&337996476\cr
524968778&-114638199&-300863016&-411132456&-492107157&-216884061&-403226667&-322669980&-332078454\cr
28705084293&19003686300&12994353171&8305795458&5569771815&3738017241&2182968630&1375264332&1161841563}\right]\)
找出因數\(p\)
for i:1 thru d do
(print("第",i,"列向量B[",i,"]=",B[ i ]),
print("產生不需要同餘p"^"rm","=p"^(r*m),"的方程式h(x)"),
printList:["h(x)=",B[ i ][1]],
for j:2 thru d do
(if B[ i ][j]>=0 then printList:append(printList,["+"]),/*若係數為正則補印+號*/
printList:append(printList,[B[ i ][j],"(",x/X,")"^(j-1)])
),
apply(print,printList),/*再用apply(print,)將全部內容印在同一行*/
print("h(x)=",hx:sum(B[ i ][j+1]*(x/X)^j,j,0,d-1),"=",factor(hx)),
print("整數解為",integerx:sublist(solve(hx,x),lambda([x],integerp(rhs(x))))),
if length(integerx)>0 then/*若有整數解*/
(for x in integerx do
(print("當整數解",x,"時，可能的因數p=P+x=",P,"+",rhs(x),"=",p: P+rhs(x)),
if mod(N,p^r)=0 then
(print("N可被p"^r,"整除，N可分解成",N,"=",p,""^r,"*",q:N/(p^r)),
i:d/*將i設為d，直接結束for迴圈*/
)
else
(print("N無法被p"^r,"整除"))
)
),
print("－－－－－－－－－")
);
第1列向量\(B[1]=\left[83412628,-91271724,-27587799,-100623168,-100709649,89252928,78155361,202680225,228782070\right]\)
產生不需要同餘\(p^{rm}=p^6\)的方程式\(h(x)\)
\(\displaystyle h(x)=83412628-91271724\left(\frac{x}{3}\right)-27587799\left(\frac{x}{3}\right)^2-100623168\left(\frac{x}{3}\right)^3-100709649\left(\frac{x}{3}\right)^4+89252928\left(\frac{x}{3}\right)^5+78155361\left(\frac{x}{3}\right)^6+202680225\left(\frac{x}{3}\right)^7+228782070\left(\frac{x}{3}\right)^8\)
\(h(x)=83412628-30423908x-3065311x^2-3726784x^3-1243329x^4+367296x^5+107209x^6+92675x^7+34870x^8\)
\(=(-2+x)^2(20853157+13247180x+7267563x^2+3024072x^3+896349x^4+232155x^5+34870x^6)\)
整數解為\([x=2]\)
當整數解\(x=2\)時，可能的因數\(p=P+x=59+2=61\)
\(N\)可被\(p^2\)整除，\(N\)可分解成\(197213=61^2 *53\)
done
