1-2.刻板訊息(Stereotyped Messages)
當大部分訊息是固定的或刻板的，假設明文m包含兩個部分。
(1)已知部分B=2kb，例如"October 19, 1995.The secret key for the day is"的ASCII值。
(2)未知部分x，例如"Squeamish Ossifrage"，x小於N31。

假設RSA加密方案採用公鑰e=3，密文C=M3=(B+x)3(modN)，若已知BCN，產生多項式p(x)=(B+x)3−C，找出x0滿足p(x0)=(B+x0)3−C0(modN)，x0N31

參考資料：
Coppersmith, D. Small Solutions to Polynomial Equations, and Low Exponent RSA Vulnerabilities. J. Cryptology 10, 233–260 (1997).
https://link.springer.com/article/10.1007/s001459900030
https://www.di.ens.fr/~fouque/ens-rennes/coppersmith.pdf

顯示數字的全部位數會造成版面凌亂，將太長的數字縮短顯示
功能表選取，編輯/設定/Worksheet/Maximum displayed number of digits調整為60



請下載LLL.zip，解壓縮後將LLL.mac放到C:\maxima-5.45.1\share\maxima\5.45.1\share目錄下
要先載入LLL.mac才能使用Coppersmith_Howgrave指令
(%i1)　load("LLL.mac");
(%o1)　C:/maxima-5.45.1/share/maxima/5.45.1/share/LLL.mac

明文訊息
(%i2)　m:"October 19,1995.The secret key for the day is";
(m)　October191995Thesecretkeyforthedayis

明文訊息轉成list
(%i3)　mlist:charlist(m);
(mlist)　[October191995Thesecretkeyforthedayis]

明文轉成數字B
(%i4)　B:0;
(B)　0

將明文訊息按照ASCII轉換成數字
(%i6)
for i:1 thru length(mlist) do
  (B:B*1000+cint(mlist[ i ])
  )$
B;
(%o6)　79099116111098101114[94digits]32100097121032105115

密碼x長度19(通常要猜測很多次，為了簡化過程設為已知)
(%i7)　length:19;
(length)　19

最後19位000為密碼x
(%i8)　B:B*1000^length;
(B)　79099116111098101114[151 digits]00000000000000000000

公鑰e
(%i9)　e:3;
(e)　3

使用RSA-230作為此次的公鑰N，因為密碼x小於上界X
https://en.wikipedia.org/wiki/RSA_numbers#RSA-230
(%i10)
N:17969491597941066732916128449573246156367561808012600070888918835531726460341490933493372247868650755230855864199929221814436684722874052065257937495694348389263171152522525654410980819170611742509702440718010364831638288518852689;
(N)　17969491597941066732[190 digits]64831638288518852689

密文C
(%i11)
C:3601065602437181695470302568875441014033597674933932563017313054187328842620219291572818766268178751411814560562534443426266396400495371407624162345901468773333852250822392143444872355448117811583996009775475129299419544905929790;
(C)　36010656024371816954[189 digits]29299419544905929790

明文M，其中x為19位密碼
(%i12)　M:B+x;
(M)　x+79099116111098101114[151 digits]00000000000000000000

產生方程式p(x)≡(B+x)^e-C(mod N)
(%i13)　px:M^e-C;
(px)　(x+79099116111098101114[151 digits]00000000000000000000)^3
　　　-36010656024371816954[189 digits]29299419544905929790

方程式p(x)同餘N
(%i14)　px:polymod(px,N);
(px)　x^3
　　　+23729734833329430334[152 digits]00000000000000000000x^2
　　　+57518597415676200644[189 digits]30736530413029557356x
　　　-78610106704918417037[189 digits]64951518554580927445

參數h
(%i15)　h:3;
(h)　3

呼叫Coppersmith_Howgrave副程式，找符合p(x)≡0(mod N)的較小的解x
執行時間需1965秒(32分鐘)
(%i18)
showtime:true$
x:Coppersmith_Howgrave(px,N,h);
showtime:false$
Evaluation took 0.0000 seconds (0.0000 elapsed) using 0 bytes.
參數h=3
p(x)最高次方k=3
\displaystyle X=ceiling(\frac{1}{\sqrt{2}} {{(hk)}^{-1/(hk-1)}} {{N}^{(h-1)/(hk-1)}} )
\displaystyle =ceiling(\frac{1}{\sqrt{2}} 9^{\frac{-1}{8}} 17969491597941066732[190 digits]64831638288518852689^{\frac{1}{4}})
=1106212689453879191977235208036134814768946283886104135857

q_{uv}=N^{h−1−v}x^u p(x)^v=
[32290262828847459190[419digits]60277520976882530721,
32290262828847459190[419digits]60277520976882530721x,
32290262828847459190[419digits]60277520976882530721x^2,
17969491597941066732[190digits]64831638288518852689(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445),
17969491597941066732[190digits]64831638288518852689x(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445),
17969491597941066732[190digits]64831638288518852689x^2(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445),
(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445)^2,
x(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445)^2,
x^2(x^3+23729734833329430334[152digits]00000000000000000000x^2+57518597415676200644[189digits]30736530413029557356x−78610106704918417037[189digits]64951518554580927445)^2]

用1106212689453879191977235208036134814768946283886104135857x取代x,得到q_{uv}=
[32290262828847459190[419digits]60277520976882530721,
35719898487071973003[476digits]68945401330960162897x,
39513804972403437655[533digits]07174654938138697729x^2,
17969491597941066732[190digits]64831638288518852689(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x^2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445),
19878079628677272620[247digits]72697393439425769473x(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x^2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445),
21989383927217453979[304digits]76931392169955293361x^2(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x^2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445),
(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445)^2,
1106212689453879191977235208036134814768946283886104135857x(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x^2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445)^2,
12237065143087845640[75digits]14878644880713124449x^2(13536796742957524728[132digits]22453052799644267793x^3+29038231098365304034[266digits]00000000000000000000x^2+63627802340810115187[246digits]11462962337597714092x−78610106704918417037[189digits]64951518554580927445)^2]

產生矩陣
M=\left[\matrix{32290262828847459190[419 digits]60277520976882530721 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\cr 0 & 35719898487071973003[476 digits]68945401330960162897 & 0 & 0 & 0 & 0 & 0 & 0 & 0\cr 0 & 0 & 39513804972403437655[533 digits]07174654938138697729 & 0 & 0 & 0 & 0 & 0 & 0\cr -14125836519472822098[419 digits]94489721632952149605 & 11433592595586423031[476 digits]00952654252887393388 & 52180224974114632454[495 digits]00000000000000000000 & 24324935533561123862[361 digits]49738464971834145377 & 0 & 0 & 0 & 0 & 0\cr 0 & -15626179606991854661[476 digits]60211676569108886485 & 12647985215283616321[533 digits]86656469202955513516 & 57722427004923821310[552 digits]00000000000000000000 & 26908552357372882659[418 digits]46542924612696483089 & 0 & 0 & 0 & 0\cr 0 & 0 & -17285878168939820519[533 digits]81089654330431192645 & 13991361741171790418[590 digits]88089882508863743212 & 63853281218922005140[609 digits]00000000000000000000 & 29766582072559977503[475 digits]49821745201359022273 & 0 & 0 & 0\cr 61795488761586594662[418 digits]61228105756354228025 & -10003576662821081972[476 digits]14874260906612109880 & 40484972307212011736[532 digits]19137568143775384464 & 36952776573071059484[552 digits]95736854560233442230 & 84321886532206985153[571 digits]51737087375795677912 & 78616926430719870707[437 digits]00000000000000000000 & 18324486606014544982[303 digits]47263486203097090849 & 0 & 0\cr 0 & 68358953819071673294[475 digits]53451685880956792425 & -11066083444337370707[533 digits]78657042708931967160 & 44784990098426820113[589 digits]33285895980663125648 & 40877630355685238072[609 digits]16675173406681041110 & 93277940880617523876[628 digits]85134090008262090584 & 86967041623524387424[494 digits]00000000000000000000 & 20270779611300936552[360 digits]36453950823767472593 & 0\cr 0 & 0 & 75619542152448801923[532 digits]51332599109148483225 & -12241441928681489685[590 digits]25369415714902456120 & 49541724343946082466[646 digits]63054918145653160336 & 45219353414264099481[666 digits]97424893201642081270 & 10318524184826785550[686 digits]19237039307576470488 & 96204045008206368850[551 digits]00000000000000000000 & 22423793631144068882[417 digits]59771418125196067201}\right]

LLL化簡
B=\left[\matrix{-35743980059088545200[416 digits]80722642160047905681 & 46922824634222774487[417 digits]96446709945328871477 & 16364540563595223391[417 digits]91327750786549963386 & -97207301925654907852[417 digits]53799084107614558498 & -68006717288291680917[417 digits]14045799605521273218 & -22083850301737135483[417 digits]51281936303530818568 & 80777391656482090011[417 digits]81364384503413938164 & -18875467394270883798[417 digits]16922840533503730786 & 0\cr -53855497698579115104[416 digits]22440422365250280904 & 74325972788911786534[417 digits]02667487696764782596 & -42360213040583351016[417 digits]52953633850159392540 & 93207789824932536615[417 digits]17899116640251965907 & 26743295137275436286[417 digits]30025365712179588699 & -19280506252185172185[417 digits]53869638507983632983 & -49522216738980310653[416 digits]09740058127199728661 & -54329652874359117291[417 digits]64740529713402819815 & 22423793631144068882[417 digits]59771418125196067201\cr 59886248145072283776[415 digits]00597316332236055074 & -13315128391293646051[417 digits]79917797014472050911 & 62558879042844919377[417 digits]09791573964335271861 & 11609344509786149688[418 digits]01042210268547527184 & -27368204832572340597[417 digits]17047835340290094121 & 17836124237718092338[417 digits]01896389102994872252 & -41032209045239051117[417 digits]16172892730225165900 & -59903924262576322974[417 digits]02209538630728871520 & -67271380893432206647[417 digits]79314254375588201603\cr 18422321777782777023[416 digits]59447344137439501243 & -21083022289790045344[417 digits]50420523383747574118 & -36876219378162584256[417 digits]55085598948970130728 & -12380026876101576575[418 digits]54404726405388353366 & 88756932911925724699[417 digits]97537235101170684352 & -14639498730317377843[418 digits]26252864402197471613 & -24797076405414559994[417 digits]60370312150796367236 & -16040594956631378629[418 digits]88876469069441734095 & 22423793631144068882[417 digits]59771418125196067201\cr -59740829062319787233[416 digits]40688909484194485472 & 69491369706263632269[417 digits]58027265444373884189 & 12960931735222433006[418 digits]05407821747361781036 & 39543384198825880926[417 digits]32633257986499319708 & 14710427912295166891[418 digits]99414848522999747063 & -49840365686044385742[417 digits]35018273347131338824 & 20876577478866948183[417 digits]89058142685532409209 & 74627815624516244194[417 digits]75977705134307053965 & -44847587262288137765[417 digits]19542836250392134402\cr -68384734181486838572[416 digits]00732056093284927437 & 91066074648706805896[417 digits]03099971016027027675 & -56710682936040096478[416 digits]32559138585047288553 & 72297457020480686195[417 digits]84446444836328742693 & -60533976691011414274[417 digits]31378847662774021892 & -12596396540470275328[418 digits]52090571729875618447 & -10880938062066902538[418 digits]67272557702486748080 & 21211756224520192677[418 digits]36650827112909148313 & 0\cr -81250564085878561445[416 digits]34678304314010262958 & 10619954225570864725[418 digits]68561477549839159629 & 29714549780697273433[417 digits]69298528953401160635 & -43469171935032168882[417 digits]37163962147628359366 & -10442288961713937897[418 digits]86086266267518395545 & 49860675747404052367[416 digits]60024965193505807977 & -13253357916405795479[418 digits]68461980318400142027 & 13768422266804784107[417 digits]71514595120672325344 & 15696655541800848217[418 digits]18399926876372470407\cr 12609739358642155334[417 digits]57802314032806547340 & -17383608959896020481[418 digits]78499216381863868700 & 81673208107847127971[417 digits]29200832824127271546 & -20653954414963469692[417 digits]04471381311357994334 & -31161557914274171485[417 digits]97962853192070185042 & -84798118358121590160[417 digits]88845627108413875920 & 41678440330915146980[417 digits]79753567469752741441 & -85952558727161375791[417 digits]14248221183260863210 & 15696655541800848217[418 digits]18399926876372470407\cr 32200663351089791530[419 digits]57114456451584344136 & 12124879742313456102[418 digits]99114197642093654073 & -25995672476988127625[417 digits]61625883063609429154 & -39995121007223712365[416 digits]35899967467362592591 & -41263422151016244631[417 digits]84020433893341684519 & -41364356553922307668[417 digits]05151574811514451551 & 75825169982584058945[417 digits]71624326376214209503 & -73205120268630001090[417 digits]81663370246906550601 & 22423793631144068882[417 digits]59771418125196067201}\right]

產生不需要同餘N^2的方程式r(x)=
-35743980059088545200[416 digits]80722642160047905681
\displaystyle +46922824634222774487[417 digits]96446709945328871477( \frac{x}{1106212689453879191977235208036134814768946283886104135857})
\displaystyle +16364540563595223391[417 digits]91327750786549963386( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^2
\displaystyle -97207301925654907852[417 digits]53799084107614558498( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^3
\displaystyle -68006717288291680917[417 digits]14045799605521273218( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^4
\displaystyle -22083850301737135483[417 digits]51281936303530818568( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^5
\displaystyle +80777391656482090011[417 digits]81364384503413938164( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^6
\displaystyle -18875467394270883798[417 digits]16922840533503730786( \frac{x}{1106212689453879191977235208036134814768946283886104135857})^7
\displaystyle +0(\frac{x}{1106212689453879191977235208036134814768946283886104135857})^8

r(x)=
-931166326910674526718530634530395511188484760951215019202x^7
+44081667002866521292[75 digits]22469878202438872436x^6
-13331579738644988581[132 digits]48018398665104447624x^5
-45414785555369089229[189 digits]85106015860347068418x^4
-71809678294997445547[246 digits]13523136087532093186x^3
+13372929188694233984[303 digits]93165215585922633914x^2
+42417543282194563777[360 digits]26324257495544510661x
-35743980059088545200[416 digits]80722642160047905681

= -(x-83113117101097109105115104032079115115105102114097103101)
(931166326910674526718530634530395511188484760951215019202x^6
-43307745643175267650[75 digits]26182397092150127034x^5
+97321380036192344236[131 digits]32115288696559115190x^4
+46223653880907932101[189 digits]89490665635112272608x^3
+75651470252841928353[246 digits]55691765552926250594x^2
-70852996826996182688[302 digits]24287941312942141920x
-43006424624419143341[360 digits]20196969789558604581)

整數解為[x=83113117101097109105115104032079115115105102114097103101]
Evaluation took 1965.2180 seconds (1977.6340 elapsed) using 132237.506 MB.}
(x)　[x=83113117101097109105115104032079115115105102114097103101]

將密碼x從個位數開始每3位數分隔出一個數字
(%i19)　makelist(mod(floor(rhs(x[1])/1000^i),1000),i,length-1,0,-1);
(%o19)　[83,113,117,101,97,109,105,115,104,32,79,115,115,105,102,114,97,103,101]

每一個數字轉換成ASCII表示法
(%i20)　makelist(ascii(i),i,%);
(%o20)　[S,q,u,e,a,m,i,s,h, ,O,s,s,i,f,r,a,g,e]

將list組合成字串，得到密碼
(%i21)　simplode(%);
(%o21)　Squeamish Ossifrage

2-1.兩個訊息有仿射關係
兩個訊息m_1和m_2=\alpha m_1+\beta有仿射關係(affine relation)
經由RSA加密，公鑰e=3，得到密文\cases{c_1=m_1^3 \pmod{N} \cr c_2=m_2^3 \pmod{N}}
若已知c_1,c_2,\alpha,\beta,N就可以回復明文
\displaystyle \frac{\beta(c_2+2\alpha^3 c_1-\beta^3)}{\alpha(c_2-\alpha^3 c_1+2\beta^3)}=\frac{3\alpha^3\beta m_1^3+3\alpha^2 \beta^2 m_1^2+3\alpha \beta^3 m_1}{3 \alpha^3 \beta m_1^2+3\alpha^2 \beta^2 m_1+3 \alpha \beta^3}=m_1 \pmod{N}
上面的式子可以從m_1^3-c_1,(\alpha m_1+\beta)^3-c_2計算歐幾里得演算法得到最大公因式。
若\alpha=\beta=1，則
\displaystyle \frac{(m+1)^3+2m^3-1}{(m+1)^3-m^3+2}=\frac{3m^3+3m^2+3m}{3m^2+3m+3}=m \pmod{N}


同樣RSA加密，換公鑰e=5，得到密文\cases{c_1=m^5 \pmod{N}\cr c_2=(m+1)^5\pmod{N}}
P(m)=c_2^3-3c_1c_2^2+3c_1^2c_2-c_1^3+37c_2^2+176c_1c_2+37c_1^2+73c_2-73c_1+14
mP(m)=2c_2^3-1c_1c_2^2-4c_1^2c_2+3c_1^3+14c_2^2-88c_1c_2-51c_1^2-9c_2+64c_1-7
\displaystyle m=\frac{mP(m)}{P(m)}
令z為未知的訊息m，則z滿足這兩個方程式\cases{z^5-c_1\equiv 0\pmod{N}\cr (z+1)^5-c_2\equiv 0\pmod{N}}
應用歐幾里得演算法找到最大公因式gcd(z^5-c_1,(z+1)^5-c_2)\in Z/N[z]，得到線性多項式z-m。

參考資料：
Coppersmith D., Franklin M., Patarin J., Reiter M. (1996) Low-Exponent RSA with Related Messages. In: Maurer U. (eds) Advances in Cryptology — EUROCRYPT ’96. EUROCRYPT 1996. Lecture Notes in Computer Science, vol 1070. Springer, Berlin, Heidelberg.
https://link.springer.com/chapter/10.1007%2F3-540-68339-9_1
https://link.springer.com/content/pdf/10.1007/3-540-68339-9_1.pdf

顯示數字的全部位數會造成版面凌亂，將太長的數字縮短顯示
功能表選取，編輯/設定/Worksheet/Maximum displayed number of digits調整為80


多項式輾轉相除法副程式
(%i1)
GCD(fx1,fx2,var):=block([temp],
fx1:expand(fx1),
fx2:expand(fx2),
while hipow(fx2,var)#1 do
  (temp:fx2,
   print(fx1,"除以",fx2,"餘式",fx2:remainder(fx1,fx2,var)),
   fx1:temp
  ),
fx2
)$

根據密文c1,c2產生方程式
(%i3)
fx1:m1^3-c1;
fx2: (alpha*m1+beta)^3-c2;
(fx1)　m_1^3-c_1
(fx2)　(\alpha m_1+\beta)^3-c_2

計算fx1和fx2輾轉相除法得到最大公因式
(%i4)　GCD:GCD(fx1,fx2,m1);
m_1^3-c_1除以\alpha^3m_1^3+3\alpha^2\beta m_1^2+3\alpha\beta^2m_1-c_2+\beta^3餘式\displaystyle -\frac{3\alpha^2\beta m_1^2+3\alpha\beta^2m_1-c_2+\alpha^3c_1+\beta^3}{\alpha^3}
\alpha^3m_1^3+3\alpha^2\beta m_1^2+3\alpha\beta^2m_1-c_2+\beta^3除以\displaystyle -\frac{3\alpha^2\beta m_1^2+3\alpha\beta^2m_1-c_2+\alpha^3c_1+\beta^3}{\alpha^3}餘式\displaystyle\frac{(\alpha c_2-\alpha^4c_1+2\alpha\beta^3)m1-\beta c_2-2\alpha^3\beta c_1+\beta^4}{3\beta}
(GCD)　\displaystyle\frac{(\alpha c_2-\alpha^4c_1+2\alpha\beta^3)m1-\beta c_2-2\alpha^3\beta c_1+\beta^4}{3\beta}

從最大公因式，解方程式得m1
(%i5)　m1:solve(GCD,m1)[1];
(m1)　\displaystyle m_1=\frac{\beta c_2+2 \alpha^3 \beta c_1-\beta^4}{\alpha c_2-\alpha^4 c_1+2 \alpha \beta^3}

範例取自Cryptanalytic Attacks on RSA
(%i10)
alpha:3;
beta:5;
N:7790302288510159542362475654705578362485767620973983941084402222135728725117099985850483876481319443405109322651368151685741199347755868542740942256445000879127232585749337061853958340278434058208881085485078737;
c1:132057584044937409231208389323398996878812486949811558724214983072091380989054308161277959733824865068687594213139826622055543700074552293693503940351187203266740911056806170880679978462212228231292575333924006;
c2:3565554769213310049242626511731772915727937147644912090997362096862084038123081043744658925329430451812652081858712220905928591327874274888835176225741122966452992998335410453929161733393892204730002674838955287;
(alpha)　3
(beta)　5
(N)　77903022885101595423624756[159 digits]78434058208881085485078737
(c1)　13205758404493740923120838[158 digits]62212228231292575333924006
(c2)　35655547692133100492426265[159 digits]93892204730002674838955287

將\alpha,\beta,c_1,c_2代入
(%i11)　m1:ev(m1,[alpha=alpha,beta=beta,c1=c1,c2=c2]);
(m1)　\displaystyle m_1=\frac{\beta c_2+2 \alpha^3 \beta c_1-\beta^4}{\alpha c_2-\alpha^4 c_1+2 \alpha \beta^3}=\frac{16978832234349095472583935[158 digits]61164325860631773696362722}{51843405275386826203986806[103 digits]08706305166524791817361975}

明文m1,m2
(%i13)
m1:ratsimp(rhs(m1)),modulus:N;
m2:alpha*m1+beta;
warning: assigning 77903022885101595423624756[159 digits]78434058208881085485078737, a non-prime, to 'modulus'
(m1)　200805001301070903002315180419000118050019172105011309190800151919090618010705
(m2)　602415003903212709006945541257000354150057516315033927572400455757271854032120

驗算
(%i15)
is(power_mod(m1,3,N)=c1);
is(power_mod(m2,3,N)=c2);
(%o14)　true
(%o15)　true

清除m1,m2,alpha,beta,N,c1,c2設定
(%i16)　kill(m1,m2,alpha,beta,N,c1,c2);
(%o16)　done

根據密文c1,c2產生方程式
(%i18)
fx1:m1^5-c1;
fx2: (m1+1)^5-c2;
(fx1)　m_1^5-c_1
(fx2)　(m_1+1)^5-c_2

計算fx1和fx2輾轉相除法得到最大公因式
(%i19)　GCD:GCD(fx1,fx2,m1)
m_1^5-c_1除以m_1^5+5m_1^4+10m_1^3+10m_1^2+5m_1-c_2+1餘式-5m_1^4-10m_1^3-10m_1^2-5m_1+c_2-c_1-1
m_1^5+5m_1^4+10m_1^3+10m_1^2+5m_1-c_2+1除以-5m_1^4-10m_1^3-10m_1^2-5m_1+c_2-c_1-1餘式\displaystyle \frac{10m_1^3+15m_1^2+(c_2-c_1+9)m_1-2c_2-3c_1+2}{5}
-5m_1^4-10m_1^3-10m_1^2-5m_1+c_2-c_1-1除以\displaystyle \frac{10m_1^3+15m_1^2+(c_2-c_1+9)m_1-2c_2-3c_1+2}{5}餘式\displaystyle \frac{(2c_2-2c_1-7)m_1^2+(-3c_2-7c_1-7)m_1+2c_2-7c_1-2}{4}
\displaystyle \frac{10m_1^3+15m_1^2+(c_2-c_1+9)m_1-2c_2-3c_1+2}{5}除以\displaystyle \frac{(2c_2-2c_1-7)m_1^2+(-3c_2-7c_1-7)m_1+2c_2-7c_1-2}{4}餘式\displaystyle \frac{4c_2^3+(148-12c_1)c_2^2+(12c_1^2+704c_1+292)c_2-4c_1^3+148c_1^2-292c_1+56)m_1-8c_2^3+(4c_1-56)c_2^2+(16c_1^2+352c_1+36)c_2-12c_1^3+204c_1^2-256c_1+28}{20c_2^2+(-40c_1-140)c_2+20c_1^2+140c_1+245}
(GCD)　\displaystyle \frac{4c_2^3+(148-12c_1)c_2^2+(12c_1^2+704c_1+292)c_2-4c_1^3+148c_1^2-292c_1+56)m_1-8c_2^3+(4c_1-56)c_2^2+(16c_1^2+352c_1+36)c_2-12c_1^3+204c_1^2-256c_1+28}{20c_2^2+(-40c_1-140)c_2+20c_1^2+140c_1+245}

從最大公因式，解方程式得m1
(%i20)　m1:solve(GCD,m1)[1];
(m1)　\displaystyle m_1=\frac{2c_2^3+(14-c_1)c_2^2+(-4c_1^2-88c_1-9)c_2+3c_1^3-51c_1^2+64c_1-7}{c_2^3+(37-3c_1)c_2^2+(3c_1^2+176c_1+73)c_2-c_1^3+37c_1^2-73c_1+14}

範例取自Cryptanalytic Attacks on RSA
(%i23)
c1:18796237015415790;
c2:7290180156009373;
N:35480779745861123;
(c1)　18796237015415790
(c2)　7290180156009373
(N)　35480779745861123

將c1,c2代入
(%i24)　m1:ev(m1,[c1=c1,c2=c2]);
(m1)　\displaystyle m_1=\frac{2c_2^3+(14-c_1)c_2^2+(-4c_1^2-88c_1-9)c_2+3c_1^3-51c_1^2+64c_1-7}{c_2^3+(37-3c_1)c_2^2+(3c_1^2+176c_1+73)c_2-c_1^3+37c_1^2-73c_1+14}=-\frac{43297460062121981374915003596042374969361963498}{7019720390981513672602639591832073102833439091}

將分數m_1同餘N
(%i25)　m1:ratsimp(rhs(m1)),modulus:N;
warning: assigning 35480779745861123, a non-prime, to 'modulus'
(m1)　-16036924398274761

明文m1,m2
(%i27)
m1:mod(m1,N);
m2:m1+1;
(m1)　19443855347586362
(m2)　19443855347586363

驗算
(%i29)
is(power_mod(m1,5,N)=c1);
is(power_mod(m2,5,N)=c2);
(%o28)　true
(%o29)　true

2-2.兩個訊息有仿射關係
上一篇文章提到兩個訊息m_1和m_2=\alpha m_1+\beta有仿射關係的話就可以回復明文，但實務上兩個訊息的\alpha和\beta值不易得知。

改成明文M向左平移k位元，各加上隨機補綴值T_1和T_2(假設T_1<T_2)後e=3次方同餘N得到密文c_1,c_2
c_1\equiv m_1^3=(2^kM+T_1)^3\pmod{N}
c_1\equiv m_2^3=(2^kM+T_2)^3\pmod{N}
若2個補綴值T_1和T_2很接近，設t為兩個補綴值的差，且t<N^{\frac{1}{9}}
t=T_2-T_1，t=(2^kM+T_2)-(2^kM+T_1)，t=m_2-m_1，m_2=m_1+t
得到兩個同餘方程式，
\cases{c_1\equiv m_1^3\pmod{N}\cr c_2\equiv m_2^3\pmod{N}}，\cases{m_1^3-c1\equiv 0\pmod{N}\cr (m_1+t)^3-c_2\equiv 0\pmod{N}}
利用Resultant消去共同變數m_1，得到t的9次方程式。
Res_{m_1}(m_1^3-c_1,(m_1+t)^3-c_2)=t^9+(3c_1-3c_2)t^6+(3c_1^2+21c_1c_2+3c_2^2)t^3+(c_1-c_2)^3\equiv 0\pmod{N}
因為t<N^{\frac{1}{9}}，可利用Coppersmith-Howgrave方法求得t較小的解。
再利用上一篇輾轉相除法求得m_1=2^kM+T_1的解，再將後面k位元刪除得到明文M。
\displaystyle m_1\equiv \frac{t(c_2+2c_1-t^3)}{c_2-c_1+2t^3}=\frac{t((m_1+t)^3+2m_1^3-t^3)}{(m_1+t)^3-m_1^3+2t^3}=\frac{t(3m_1^3+3m_1^2t+3m_1t^2)}{3m_1^2t+3m_1t^2+3t^3}\pmod{N}

本篇文章先介紹要如何計算Resultant(結式)。

參考資料：
D. Coppersmith. Finding a small root of a univariate modular equation. In Proceedings of Eurocrypt 1996, volume 1070 of Lecture Notes in Computer Science, pages 155–165. Springer, 1996.
https://link.springer.com/chapter/10.1007/3-540-68339-9_14
https://isc.tamu.edu/resources/preprints/1996/1996-02.pdf
https://canvas.mit.edu/courses/7 ... load?download_frd=1

公式	範例
設f(x)=a_nx^n+a_{n-1}x^{n-1}+\ldots+a_0和 g(x)=b_mx^m+b_{m-1}x^{m-1}+\ldots+b_0 是次數為n和m的多項式 設f(x)和g(x)的Sylvester矩陣為 Syl(f,g)=\left[\matrix{a_n&a_{n-1}&\ldots&\ldots&a_0&0&0\cr 0&a_n&a_{n-1}&\ldots&a_1&a_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&a_n&a_{n-1}&\ldots&a_0\cr b_m&b_{m-1}&\ldots&\ldots&b_0&0&0\cr 0&b_m&b_{m-1}&\ldots&b_1&b_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&b_m&b_{m-1}&\ldots&b_0}\right] \matrix{　\cr 　 \cr f(x)取m列\cr 　\cr 　\cr　\cr g(x)取n列\cr　} 則f(x)和g(x)的Resultant為Sylvester矩陣的行列式值 Res(f,g)=det(Syl(f,g)) Res(f,g)=\left|\matrix{a_n&a_{n-1}&\ldots&\ldots&a_0&0&0\cr 0&a_n&a_{n-1}&\ldots&a_1&a_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&a_n&a_{n-1}&\ldots&a_0\cr b_m&b_{m-1}&\ldots&\ldots&b_0&0&0\cr 0&b_m&b_{m-1}&\ldots&b_1&b_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&b_m&b_{m-1}&\ldots&b_0}\right|	f(x)=1x^3+2x^2+3x+4 g(x)=5x^2+6x+7 f(x)和g(x)的Sylvester矩陣為 Syl(f,g)=\left[\ \matrix{1&2 &3 &4 &0\cr 0 &1&2 &3 &4 \cr 5&6&7&0&0 \cr 0&5&6&7&0 \cr 0&0&5&6&7} \right] \matrix{f(x)取2列\cr　\cr 　 \cr g(x)取3列 \cr　} f(x)和g(x)的Resultant為 Res(f,g)= \left|\ \matrix{1&2 &3 &4 &0\cr 0 &1&2 &3 &4 \cr 5&6&7&0&0 \cr 0&5&6&7&0 \cr 0&0&5&6&7} \right|=832 超過3階行列式不能用交叉相乘方式計算，而降階方式計算量又太大，實務上可採用PA=LU分解，再個別計算矩陣P,L,U的行列式值，det(P)\cdot det(A)=det(L)\cdot det(U)，進而得到矩陣A的行列式值。 參考資料：https://ccjou.wordpress.com/2012/04/13/palu-%E5%88%86%E8%A7%A3/ P=\left[ \matrix{1&0&0&0&0\cr 0&1&0&0&0\cr 0&0&0&1&0\cr 0&0&1&0&0\cr 0&0&0&0&1}\right]，det(P)=-1 L=\left[ \matrix{\displaystyle 1&0&0&0&0\cr 0&1&0&0&0\cr 0&5&1&0&0\cr 5&-4&0&1&0\cr 0&0&-\frac{5}{4}&\frac{1}{2}&1}\right]，det(L)=1\cdot 1\cdot 1\cdot 1\cdot 1=1 U=\left[ \matrix{1&2&3&4&0\cr 0&1&2&3&4\cr 0&0&-4&-8&-20\cr 0&0&0&-8&16\cr 0&0&0&0&-26}\right]，det(U)=1\cdot 1\cdot (-4)\cdot(-8)\cdot(-26)=-832 \displaystyle Res(f,g)=\frac{det(L)\cdot det(U)}{det(P)}=\frac{1\cdot (-832)}{-1}=832
判斷是否有共同根 若Res(f,g)=0，則f,g有共同根 [證明] 設x=x_0是f(x)=0,g(x)=0的共同根，則f(x_0)=g(x_0)=0 \left[\matrix{a_n&a_{n-1}&\ldots&\ldots&a_0&0&0\cr 0&a_n&a_{n-1}&\ldots&a_1&a_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&a_n&a_{n-1}&\ldots&a_0\cr b_m&b_{m-1}&\ldots&\ldots&b_0&0&0\cr 0&b_m&b_{m-1}&\ldots&b_1&b_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&b_m&b_{m-1}&\ldots&b_0}\right] \left[\matrix{x_0^{n+m-1}\cr x_0^{n+m-2}\cr \vdots \cr x_0^n \cr x_0^{n-1}\cr \vdots \cr x_0 \cr 1}\right]=\left[\matrix{f(x_0)x_0^{m-1}\cr f(x_0)x_0^{m-2}\cr \vdots \cr f(x_0)x_0 \cr f(x_0)\cr g(x_0)x_0^{n-1} \cr g(x_0)x_0^{n-2}\cr \vdots \cr g(x_0)x_0 \cr g(x_0)} \right]=0 Syl(f,g)x=0是齊次方程組 若x有非0解，則Syl(f,g)的行列式值要為0，得到 Res(f,g)=\left|\matrix{a_n&a_{n-1}&\ldots&\ldots&a_0&0&0\cr 0&a_n&a_{n-1}&\ldots&a_1&a_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&a_n&a_{n-1}&\ldots&a_0\cr b_m&b_{m-1}&\ldots&\ldots&b_0&0&0\cr 0&b_m&b_{m-1}&\ldots&b_1&b_0&0\cr 0&0&\ddots&\ddots&\ddots&\ddots&\vdots\cr 0&0&0&b_m&b_{m-1}&\ldots&b_0}\right|=0	f(x)=(x^2-1)(x^2+x+2)=x^4+x^3+x^2-x-2 g(x)=(x^2-1)(x^3+x^2+x-1)=x^5+x^4-2x^2-x+1有共同根 設x=x_0是f(x)=0,g(x)=0的共同根，則f(x_0)=g(x_0)=0 \left[\matrix{1&1&1&-1&-2&0&0&0&0\cr 0&1&1&1&-1&-2&0&0&0\cr 0&0&1&1&1&-1&-2&0&0\cr 0&0&0&1&1&1&-1&-2&0\cr 0&0&0&0&1&1&1&-1&-2\cr 1&1&0&-2&-1&1&0&0&0\cr 0&1&1&0&-2&-1&1&0&0\cr 0&0&1&1&0&-2&-1&1&0\cr 0&0&0&1&1&0&-2&-1&1}\right] \left[\matrix{x_0^8\cr x_0^7\cr x_0^6\cr x_0^5\cr x_0^4\cr x_0^3\cr x_0^2\cr x_0^1\cr 1} \right]= \left[\matrix{f(x_0)\cdot x_0^4\cr f(x_0)\cdot x_0^3\cr f(x_0)\cdot x_0^2\cr f(x_0)\cdot x_0\cr f(x_0)\cr g(x_0)\cdot x_0^3\cr g(x_0)\cdot x_0^2\cr g(x_0)\cdot x_0\cr g(x_0)} \right]=0 Syl(f,g)x=0是齊次方程組 若x有非0解，則Syl(f,g)的行列式值要為0，得到 Res(f,g)=\left|\ \matrix{1&1&1&-1&-2&0&0&0&0\cr 0&1&1&1&-1&-2&0&0&0\cr 0&0&1&1&1&-1&-2&0&0\cr 0&0&0&1&1&1&-1&-2&0\cr 0&0&0&0&1&1&1&-1&-2\cr 1&1&0&-2&-1&1&0&0&0\cr 0&1&1&0&-2&-1&1&0&0\cr 0&0&1&1&0&-2&-1&1&0\cr 0&0&0&1&1&0&-2&-1&1}\right|=0
兩變數方項式消除共同變數 f(x,y)=0,g(x,y)=0為兩變數方程式，已知f(x)和g(x)有共同根，可利用Resultant消除共同變數	設f(x,y)=x^2y^2-25x^2+9=0 g(x,y)=4x+y=0 若要消去x變數，先以x為變數降冪排列 f(x,y)=(y^2-25)x^2+0x+9，g(x,y)=4x+y 計算f(x,y)和g(x,y)的Resultant Res_x(f,g)=\left|\ \matrix{y^2-25&0&9 \cr 4&y&0 \cr 0&4&y} \right|\ =y^4-25y^2+144=0 解得y=\pm 3,\pm 4，代回g(x,y)=0 得\displaystyle (x,y)=(1,-4),(-1,4),(\frac{3}{4},-3),(-\frac{3}{4},3) 若要消去y變數，先以y為變數降冪排列 f(x,y)=x^2y^2+0y+(-25x^2+9),g(x,y)=y+4x 計算f(x,y)和g(x,y)的Resultant Res_y(f,g)=\left|\ \matrix{x^2&0&-25x^2+9 \cr 1&4x&0 \cr 0&1&4x} \right|\ =16x^4-25x^2+9=0 解得\displaystyle x=\pm \frac{3}{4},\pm 1，代回g(x,y)=0 得\displaystyle (x,y)=(1,-4),(-1,4),(\frac{3}{4},-3),(-\frac{3}{4},3) 範例出處 http://buzzard.ups.edu/courses/2 ... ts-ups-434-2016.pdf

取多項式P的x的各項係數
(%i1)
coeffs(P,x):=block (local(l),l:[],
  for i from hipow(P,x) step -1 thru 0 do l:cons(coeff(P,x,i),l),
  reverse(l));
(%o1)　coeffs(P,x):=block(local(l),l:[],for i from hipow(P,x) step -1 thru 0 do l:cons(coeff(P,x,i),l),reverse(l))

計算多項式P和Q的Sylvester矩陣
副程式出處http://www.lpthe.jussieu.fr/~talon/subresultant.mac
(%i2)
result(P,Q,x):=block(local(mat,len1,len2,ll1,ll2),
  len1:hipow(P,x)+1,
  len2:hipow(Q,x)+1,   /* assume len1 >= len2 */
  mat:zeromatrix(len1+len2-2,len1+len2-2),
  ll1:coeffs(P,x),
  ll2:coeffs(Q,x),
  for i from 1 thru len2-1 do (
    for j from i thru i+len1-1 do (
      mat[i,j]:ll1[j-i+1])),
  for i from len2 thru len2+len1-2 do (
    for j from i-len2+1 thru i do (
       mat[i,j]:ll2[j-i+len2])),
  mat)$

第2個例子
(%i5)
f:x^3+2*x^2+3*x+4;
g:5*x^2+6*x+7;
result(f,g,x);
(f)　x^3+2x^2+3x+4
(g)　5x^2+6x+7
(%o5)　\left[\matrix{1&2&3&4&0\cr 0&1&2&3&4\cr 5&6&7&0&0\cr 0&5&6&7&0\cr 0&0&5&6&7}\right]

計算LU分解
(%i6)　[P,L,U]:get_lu_factors(lu_factor(%,generalring));
(%o6)　[\left[\matrix{1&0&0&0&0\cr 0&1&0&0&0\cr 0&0&0&1&0\cr 0&0&1&0&0\cr 0&0&0&0&1} \right],\left[\matrix{\displaystyle 1&0&0&0&0\cr 0&1&0&0&0\cr 0&5&1&0&0\cr 5&-4&0&1&0\cr 0&0&-\frac{5}{4}&\frac{1}{2}&1} \right],\left[\matrix{1&2&3&4&0\cr 0&1&2&3&4\cr 0&0&-4&-8&-20\cr 0&0&0&-8&16\cr 0&0&0&0&-26} \right]]

計算矩陣P,L,U的行列式值
(%i9)
detP:determinant(P);
detL:determinant(L);
detU:determinant(U);
(detP)　-1
(detL)　1
(detU)　-832

計算行列式值
(%i10)　detL*detU/detP;
(%o10)　832

第3個例子
(%i14)
f:x^4+x^3+x^2-x-2;
g:x^5+x^4-2x^2-x+1;
result(f,g,x);
determinant(%);
(f)　x^4+x^3+x^2-x-2
(g)　x^5+x^4-2*x^2-x+1
(%o13)　\left[\matrix{1&1&1&-1&-2&0&0&0&0\cr 0&1&1&1&-1&-2&0&0&0\cr 0&0&1&1&1&-1&-2&0&0\cr 0&0&0&1&1&1&-1&-2&0\cr 0&0&0&0&1&1&1&-1&-2\cr 1&1&0&-2&-1&1&0&0&0\cr 0&1&1&0&-2&-1&1&0&0\cr 0&0&1&1&0&-2&-1&1&0\cr 0&0&0&1&1&0&-2&-1&1}\right]
(%o14)　0

第4個例子
(%i16)
f:x^2*y^2-25*x^2+9;
g:4*x+y;
(f)　x^2y^2-25x^2+9
(g)　y+4x

消掉共同變數x
(%i19)
result(f,g,x);
determinant(%);
ratsimp(%);
(%o17)　\left[\matrix{y^2-25&0&9\cr 4&y&0\cr 0&4&y}\right]
(%o18)　y^2(y^2-25)+144
(%o19)　y^4-25y^2+144

先解出y
(%i20)　solve(%,y);
(%o20)　[y=-4,y=4,y=-3,y=3]

將y代回g(x,y)=0求x
(%i21)　map(lambda([y],[solve(ev(g,y),x)[1],y]),%);
(%o21)　\displaystyle [[x=1,y=-4],[x=-1,y=4],[x=\frac{3}{4},y=-3],[x=-\frac{3}{4},y=3]]

消掉共同變數y
(%i23)
result(f,g,y);
determinant(%);
(%o22)　\left[\matrix{x^2&0&9-25x^2\cr 1&4x&0\cr 0&1&4x}\right]
(%o23)　16x^4-25x^2+9

先解出x
(%i24)　solve(%,x);
(%o24)　　\displaystyle [x=-1,x=1,x=-\frac{3}{4},x=\frac{3}{4}]

將x代回g(x,y)=0求y
(%i25)　map(lambda([x],[x,solve(ev(g,x),y)[1]]),%);
(%o25)　\displaystyle [[x=-1,y=4],[x=1,y=-4],[x=-\frac{3}{4},y=3],[x=\frac{3}{4},y=-3]]

第1個例子
(%i30)
f:m1^3-c1;
g: (m1+t)^3-c2;
result(f,expand(g),m1);
determinant(%);
ratsimp(%);
(f)　m1^3-c1
(g)　(t+m1)^3-c2
(%o28)　\left[\matrix{1&0&0&-c1&0&0\cr 0&1&0&0&-c1&0\cr 0&0&1&0&0&-c1\cr 1&3t&3t^2&t^3-c2&0&0\cr 0&1&3t&3t^2&t^3-c2&0\cr 0&0&1&3t&3t^2&t^3-c2}\right]
(%o29)　(t^3-c2)^3+c1((t^3-c2)^2+c1(10t^3-c2)-9t^3(t^3-c2))+c1((t^3-c2)^2-c1(8t^3+c2)+c1(t^3-c2+c1))+c1(3t^2(9t^4-3t(t^3-c2))-(t^3-c2)(8t^3+c2))
(%o30)　t^9+(3c1-3c2)t^6+(3c2^2+21c1c2+3c1^2)t^3-c2^3+3c1c2^2-3c1^2c2+c1^3

也可以使用maxima內建指令計算resultant
(%i31)　resultant(f,g,m1);
(%o31)　t^9+c1(3t^6+21c2t^3+3c2^2)-3c2t^6+c1^2(3t^3-3c2)+3c2^2t^3-c2^3+c1^3

整理成t的多項式
(%i32)　ratsimp(%);
(%o32)　t^9+(3c1-3c2)t^6+(3c2^2+21c1c2+3c1^2)t^3-c2^3+3c1c2^2-3c1^2c2+c1^3

2-3.兩個訊息有仿射關係
設明文M乘10^5倍，各加上隨機補綴值T_1和T_2(假設0\le T_1,T_2<10^5,T_1<T_2)後e=3次方同餘N得到密文c_1,c_2。
c1=1881676371789154860897069000
c2=1881678004162711039676405223
N=54957464841358314276864542898551
若2個補綴值T_1和T_2很接近，設t為兩個補綴值的差，t未知但足夠小(t<N^{\frac{1}{9}})。
將密文c_1,c_2，公鑰N代入9次同餘方程式
t^9+(3c_1-3c_2)t^6+(3c_1^2+21c_1c_2+3c_2^2)t^3+(c_1-c_2)^3\equiv 0\pmod{N}

雖然論文說能找到比N^{\frac{1}{9}}=3362還小的解，但以Coppersmith_Howgrave()副程式精算的上界X卻小很多。
\displaystyle X=ceiling(\frac{1}{\sqrt{2}}(hk)^{-1/(hk-1)}N^{(h-1)/(hk-1)})
設h=3，9次方程式k=9，計算出X=172
而且會產生27 \times 27的超大矩陣，導致LLL化簡的時間大幅增加還無法求得較小的解。

換另一個作法，令t^3=x，先解3次同餘方程式，不僅Coppersmith_Howgrave()副程式精算的上界X提高很多。
x^3+(3c_1-3c_2)x^2+(3c_1^2+21c_1c_2+3c_2^2)x+(c_1-c_2)^3\equiv 0\pmod{N}
設h=3，3次方程式k=3，計算出X=46260610
而且產生9\times 9較小的矩陣，LLL化簡的時間變得更短了

所得到的較小x的解，再解一次3次同餘方程式，得到兩個補綴值的差t。
t^3-x\equiv 0\pmod{N}

再將密文c_1,c_2，補綴值的差t，公鑰N代入以下公式得到有補綴值的m_1，再將後5位補綴值刪除得到明文M。
\displaystyle m_1\equiv \frac{t(c_2+2c_1-t^3)}{c_2-c_1+2t^3}\pmod{N}


請下載LLL.zip，解壓縮後將LLL.mac放到C:\maxima-5.46.0\share\maxima\5.46.0\share目錄下
要先載入LLL.mac才能使用Coppersmith_Howgrave指令
(%i1)　load("LLL.mac");
(%o1)　C:/maxima-5.46.0/share/maxima/5.46.0/share/LLL.mac

已知密文c_1,c_2，公鑰N
(%i4)
c1:1881676371789154860897069000;
c2:1881678004162711039676405223;
N:54957464841358314276864542898551;
(c1)　1881676371789154860897069000
(c2)　1881678004162711039676405223
(N)　54957464841358314276864542898551

9次同餘方程式
(%i5)　fx:x^9+(3*c1-3*c2)*x^6+(3*c1^2+21*c1*c2+3*c2^2)*x^3+(c1-c2)^3;
(fx)　x^9-4897120668536338008669x^6+95599144073213399280057649148373498671072396308207166187x^3
-4349693466736349905369332025355097289037604766707550608234921567

同餘方程式係數同餘N，讓係數變小
(%i6)　fx:polymod(fx,N);
(fx)　x^9-4897120668536338008669x^6-1516451737447758219766669752498x^3-21000738238808374545647388458802

參數h
(%i7)　h:3;
(h)　3

利用Coppersmith_Howgrave方法解9次同餘方程式，執行時間需2262秒(37分鐘)
但執行結果無法求得較小的整數解，故省略執行過程
(%i10)
showtime:true$
Coppersmith_Howgrave(fx,N,h);
showtime:false$
...
整數解為[]
Evaluation took 2262.6870 seconds (2264.0150 elapsed) using 420496.373 MB.
(%o9)　[]

令t^3=x，化簡成3次同餘方程式
(%i11)　fx:x^3+(3*c1-3*c2)*x^2+(3*c1^2+21*c1*c2+3*c2^2)*x+(c1-c2)^3;
(fx)　x^3-4897120668536338008669x^2+95599144073213399280057649148373498671072396308207166187x
-4349693466736349905369332025355097289037604766707550608234921567

同餘方程式係數同餘N，讓係數變小
(%i12)　fx:polymod(fx,N);
(fx)　x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802

利用Coppersmith_Howgrave方法解3次同餘方程式，執行時間需43秒
(%i15)
showtime:true$
X:Coppersmith_Howgrave(fx,N,h);
showtime:false$
Evaluation took 0.0000 seconds (0.0000 elapsed) using 0 bytes.
參數h=3
p(x)最高次方k=3
\displaystyle X=ceiling(\frac{1}{\sqrt{2}}(hk)^{-1/(hk-1)}N^{(h-1)/(hk-1)})=ceiling(\frac{1}{\sqrt{2}}9^{-1/8}54957464841358314276864542898551^{1/4})=46260610
q_uv=N^{h-1-v} x^u p(x)^v = [3020322941789135243826751301254310584993059920451586964677899601,
3020322941789135243826751301254310584993059920451586964677899601x,
3020322941789135243826751301254310584993059920451586964677899601x^2,
54957464841358314276864542898551(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802),
54957464841358314276864542898551x(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802),
54957464841358314276864542898551x^2(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802),
(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802)^2,
x(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802)^2,
x^2(x^3-4897120668536338008669x^2-1516451737447758219766669752498x-21000738238808374545647388458802)^2]
用46260610x取代x,得到q_{uv}= [3020322941789135243826751301254310584993059920451586964677899601,
139721981684159887751924249514318172791235797686641888454048089061016610x,
6463624103118063744935544456324562407407970654820642611436221569296955598732100x^2,
54957464841358314276864542898551(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802),
2542365847614788847019462641858137376110x(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802),
117611394953827177084317023684568968482648027100x^2(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802),
(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802)^(2),
46260610x(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802)^(2),
2140044037572100x^2(98999742604948264981000x^3-10480053887972286407738186731512534900x^2-70151982409893138378920200419106503780x-21000738238808374545647388458802)^(2)]
產生矩陣M=\left[ \matrix{ 3020322941789135243826751301254310584993059920451586964677899601&0&0&0&0&0&0&0&0\cr 0&139721981684159887751924249514318172791235797686641888454048089061016610&0&0&0&0&0&0&0\cr 0&0&6463624103118063744935544456324562407407970654820642611436221569296955598732100&0&0&0&0&0&0\cr -1154147333401880370164435877744830779893335235336308145728995902&-3855375106843289059874469643405077531431008453876682145059536838022780&-575957193083777435966023818707724430220694882620408968134835546929900&5440774873514967046498526949032608069078072148942531000&0&0&0&0&0\cr 0&-53391559673044361070832604010361296224641422921151169969392245114020220&-178352004221385726316119489130401363701292623991492180806562660444404996695800&-26644131085943325292024201127948743853911779893898517315388054650660801239000&251693564521475219700920220763687359166473755234092339003910000&0&0&0&0\cr 0&0&-2469926119326432700196969469407759883742609255600435024997766588244094929534200&-8250672510003878744676740400060736629653654574347063094341880675381046234195692438000&-1232573756955700653437467678941596939415709824077480689105413794852905568404895790000&11643497827837801763248586973862843084330167466119804398647668985100000&0&0&0\cr 441031006574948269107459443949256719868815747806770028851275204&2946486839047310595109751478204161934049771268347937188112057574543120&4921301076215693158799367543245210199690179969531456795727185981205628668000&1470393112007528062187159841915196407503917007629683675986668119372838320000&109831529494789146606929429768530109268611975609482301299443587987461650000&-2075045274790487347904810593080778841186998131216420673800000&9800949035846008678895673107522190930361000000&0&0\cr 0&20402363393071117647355229427353424907730536471147343664377590214914440&136306278531300406989240120330126235607922189234249266563748531753485202503200&227662389779394457098885610164724803415889536400306605558985017074020917615167480000&68021282301266572748895948454500554080939778162323820958185619069740318114575200000&5080873551661937743415985648044361658132643844999773042316053040888648280606500000&-95992860189425566911358759970378608468403657618931662386599018000000&453397880977148227551007964314572140974967380210000000&0\cr 0&0&943825776005139675748417799999320121820808332782523517793742593671973092208400&6305611591687860920570511402945041036226201306511803963291610965520595093771558952000&10531801025252553016013278646442369888149133646495387760188037870704622801637392876962800000&3146706012238795427973303402033752877122263511053778975056451231393819657574298642872000000&235044309832747753792446979829777477365827565182434950799096426463863811536307859965000000&-4440688268007542274915272165073296358699518927682926250318126398080980000000&20974462546710273066928434544050339130507985738616528100000000}\right]
LLL化簡B=\left[ \matrix{ 53771676921999198037128719975410481845009948321580737058118762&50209380384868519769271355772363021078882998018525986723096380&107113429339000465042631199484653568414097906566121374434204800&2566158289017059510580947933126354061025546078929008772549000&-39201258910319328503256252507864134476510009966583646884840000&-82244768124108488515911014612606850013994967332618611133000000&-103326720791445976823937028600478023227439055205294880233000000&37553747697339884108136344819118362708215810289041024760000000&-41948925093420546133856869088100678261015971477233056200000000\cr 145771311156080366680612542702670552709050420859148668552121201&-102701813735090026820390486307494510321459974199722868939788630&55690597658188982675500688460805070018767816245101070397610800&-3690675191418901829302352712732398605405813417421315444807000&-34161730786036250811133667107244245173497899683329768637600000&105158149470961723712865498187081390865505368625840232003900000&-83943991120073620894214705274131521444270306664078709465000000&-94290002110963796909267598874077903051034921264274818610000000&0\cr 22520257569452314079351935267514264370266330796992552972684997&-98310513164562058062910678303618723610419189844159039543246600&-59216755675973590492737933555662165421890387526607433485560200&210122354709341487106669358644005820317705567891681532969573000&-13547462576968605065539206672838828517502007656710102017580000&23570869821207231876579884845247791136173699754579501358900000&-40641684543556222622266012007342445503983888770262907798000000&-86958932882212997194203757223689795368058144049335172070000000&41948925093420546133856869088100678261015971477233056200000000\cr 169015127177354263666121860022442922657381338158956682760146316&-44086339219025359808922118654836582121338928607019431900103130&-55743427083965516899322695992901696889154913782390199199450100&-99213756448548832412671199081590293096668700399267187647555000&-13681783036472695133702522618982681148717527439219162384970000&-15336537742431604945023054843850719026627918786732999429600000&64672345006900325309865538762091565046333246270689577751000000&-71951332896677550420724086465591699219582263948746051800000000&62923387640130819200785303632151017391523957215849584300000000\cr -111109391058117478677425543043666767886469428988351808979162264&68706199548070844584789292006960471026818898889658819004932390&114242801459017508543027293838247104016630443508968997031447000&-70229925147280191228538106815486527561582710206602914483698000&-108177296750967223863314840883237590333571940265343423758290000&21114534120011920453814163485385935919207972703328940545600000&-34780235497393011659779194620008032778904564753918395174000000&44536869417475438776116536835516252809621296087125152680000000&83897850186841092267713738176201356522031942954466112400000000\cr 32124330345979841767629832716536687416322338071176926036776624&27164161165425888332075134727974947827481095848187898522695730&-79639511287230298148663766719580830450559479484044659434682900&-51041393060884297450399678845972688262211525661235763838740000&-120857089745132822141399493042839109949450984516626943707450000&220466029767922832433728342036436413650585457627045150329700000&-42081776027182400361358030886271577560506458471894525675000000&58148929768074224345419762221788597120206147540057861770000000&-41948925093420546133856869088100678261015971477233056200000000\cr -118279971490860306266688062364858934470851326591670377750842105&186394623214470748860705331927751369982342478110034345905919190&-5844542338433954468770696410970613047578620794759879329291700&-31690810413815853326709322936928777204830797244545090618294000&36631940618043663684819094023611438217838817287345595588310000&14670122074875238001436616319545787277339265229080644956400000&-87941416624880050595555558070529584290459598827529416836000000&-58985200396438287348022673963584538334771708074420956690000000&62923387640130819200785303632151017391523957215849584300000000\cr -67402575181872320450341662740231249645051561595949203267092646&81687875727598274014317606405009186128237310267665643776123970&-60606071715445114152730683837370692931748434259046928181917500&118448345919188327186093480437680808892299725583302023014092000&5659950991122719764062671142382391292964286805841933824120000&61863373560034972019278274616457949997753618484681115789000000&87493980444210846078712251695585839215963143393159467086000000&-110345238655812995600366365368792036283900354142468092020000000&-125846775280261638401570607264302034783047914431699168600000000\cr 391536656683468813127560730295755766606487863741754115343991921&447549940911495565399271813369673072505177321547345359598698170&345431307115038022995995779076001023195500573892855741305309800&330776015310331719481270882134135987680304349298242517737362000&334129783160230133456909028931418304394164386380225530980110000&298169048606388347296126093056185837124427112397441890970200000&317606145581069503748397520803157118676806368270962392206000000&313838224166290677083070069124989059103490292778965585030000000&440463713480915734405497125425057121740667700510947090100000000}\right]
產生不需要同餘N^2的方程式
r(x)=53771676921999198037128719975410481845009948321580737058118762
\displaystyle +50209380384868519769271355772363021078882998018525986723096380 ( \frac{x}{46260610} )
\displaystyle +107113429339000465042631199484653568414097906566121374434204800 ( \frac{x}{46260610} )^2
\displaystyle +2566158289017059510580947933126354061025546078929008772549000 (\frac{x}{46260610})^3
\displaystyle -39201258910319328503256252507864134476510009966583646884840000 (\frac{x}{46260610})^4
\displaystyle -82244768124108488515911014612606850013994967332618611133000000 (\frac{x}{46260610})^5
\displaystyle -103326720791445976823937028600478023227439055205294880233000000 (\frac{x}{46260610})^6
\displaystyle +37553747697339884108136344819118362708215810289041024760000000 (\frac{x}{46260610})^7
\displaystyle -41948925093420546133856869088100678261015971477233056200000000 (\frac{x}{46260610})^8

r(x)= -2x^8+82827356x^7-10542521995935153x^6-388196401064276818832330x^5-8559622143684325599618510860324x^4
+25920858191087824461843813370140390129x^3+50051974379238314010577044980314319102425126688x^2+1085359237261863165428889843267588150672526756965072158x+53771676921999198037128719975410481845009948321580737058118762
= -(x-45499293)(2x^7+8171230x^6+10914307183875543x^5+884789661515435025323429x^4+48816926196345927839271626696021x^3
+2195214770175831075654032836258741023024x^2+49828745646919485619110961668243162515263595344x+1181813460749801058164326201187576295856068539825019634)
整數解為[x=45499293]
Evaluation took 43.1100 seconds (43.1220 elapsed) using 5541.766 MB.
(X)　[x=45499293]

再解一次3次同餘方程式
(%i16)　fx:x^3-rhs(X[1]);
(fx)　x^3-45499293

利用Coppersmith_Howgrave方法解3次同餘方程式，執行時間需0.328秒
(%i19)
showtime:true$
X:Coppersmith_Howgrave(fx,N,h);
showtime:false$
Evaluation took 0.0000 seconds (0.0000 elapsed) using 0 bytes.
參數h=3
p(x)最高次方k=3
\displaystyle X=ceiling(\frac{1}{\sqrt{2}}(hk)^{-1/(hk-1)}N^{(h-1)/(hk-1)})=ceiling(\frac{1}{\sqrt{2}}9^{-1/8}54957464841358314276864542898551^{1/4})=46260610
q_{uv}=N^{h-1-v} x^u p(x)^v =
[3020322941789135243826751301254310584993059920451586964677899601,
3020322941789135243826751301254310584993059920451586964677899601x,
3020322941789135243826751301254310584993059920451586964677899601x^2,
54957464841358314276864542898551(x^3-45499293),
54957464841358314276864542898551x(x^3-45499293),
54957464841358314276864542898551x^2(x^3-45499293),
(x^3-45499293)^2,
x(x^3-45499293)^2,
x^2(x^3-45499293)^2]
用46260610x取代x,得到q_uv=
[3020322941789135243826751301254310584993059920451586964677899601,
139721981684159887751924249514318172791235797686641888454048089061016610x,
6463624103118063744935544456324562407407970654820642611436221569296955598732100x^2,
54957464841358314276864542898551(98999742604948264981000x^3-45499293),
2542365847614788847019462641858137376110x(98999742604948264981000x^3-45499293),
117611394953827177084317023684568968482648027100x^2(98999742604948264981000x^3-45499293),
(98999742604948264981000x^3-45499293)^2,
46260610x(98999742604948264981000x^3-45499293)^2,
2140044037572100x^2(98999742604948264981000x^3-45499293)^2]
產生矩陣M=\left[ \matrix{ 3020322941789135243826751301254310584993059920451586964677899601&0&0&0&0&0&0&0&0\cr 0&139721981684159887751924249514318172791235797686641888454048089061016610&0&0&0&0&0&0&0\cr 0&0&6463624103118063744935544456324562407407970654820642611436221569296955598732100&0&0&0&0&0&0\cr -2500525795354160459269142958652241224443&0&0&5440774873514967046498526949032608069078072148942531000&0&0&0&0&0\cr 0&-115675848613818628883670707444457456909880090230&0&0&251693564521475219700920220763687359166473755234092339003910000&0&0&0&0\cr 0&0&-5351235319142904201522225965512143075699768000894840300&0&0&11643497827837801763248586973862843084330167466119804398647668985100000&0&0&0\cr 2070185663499849&0&0&-9008836591414248716424316866000&0&0&9800949035846008678895673107522190930361000000&0&0\cr 0&95768051606757749647890&0&0&-416754276109143908313505917054448260000&0&0&453397880977148227551007964314572140974967380210000000&0\cr 0&0&4430288485840093620938676612900&0&0&-19279307032917423776366854961548179721038600000&0&0&20974462546710273066928434544050339130507985738616528100000000}\right]
LLL化簡B=\left[ \matrix{ 2070185663499849&0&0&-9008836591414248716424316866000&0&0&9800949035846008678895673107522190930361000000&0&0\cr 0&95768051606757749647890&0&0&-416754276109143908313505917054448260000&0&0&453397880977148227551007964314572140974967380210000000&0\cr -2500525795354160459269142958652241224443&0&0&5440774873514967046498526949032608069078072148942531000&0&0&0&0&0\cr 0&0&4430288485840093620938676612900&0&0&-19279307032917423776366854961548179721038600000&0&0&20974462546710273066928434544050339130507985738616528100000000\cr 0&-115675848613818628883670707444457456909880090230&0&0&251693564521475219700920220763687359166473755234092339003910000&0&0&0&0\cr 3020322941789135243826751301254310584993059920451586964677899601&0&0&0&0&0&0&0&0\cr 0&0&-5351235319142904201522225965512143075699768000894840300&0&0&11643497827837801763248586973862843084330167466119804398647668985100000&0&0&0\cr 0&139721981684159887751924249514318172791235797686641888454048089061016610&0&0&0&0&0&0&0\cr 0&0&6463624103118063744935544456324562407407970654820642611436221569296955598732100&0&0&0&0&0&0}\right]
產生不需要同餘N^2的方程式
\displaystyle r(x)= 2070185663499849 + 0 ( \frac{x}{46260610} ) + 0 (\frac{x}{46260610})^2
\displaystyle -9008836591414248716424316866000 (\frac{x}{46260610})^3 + 0 (\frac{x}{46260610})^4 + 0 (\frac{x}{46260610})^5
\displaystyle +9800949035846008678895673107522190930361000000 (\frac{x}{46260610})^6 + 0 (\frac{x}{46260610})^7 + 0 (\frac{x}{46260610})^8
r(x)=x^6-90998586x^3+2070185663499849=(x-357)^2(x^2+357x+127449)^2
整數解為 [x=357]
Evaluation took 0.3280 seconds (0.3280 elapsed) using 88.429 MB.
(X)　[x=357]

從密文c_1,c_2，兩個補綴值的差t=T_2-T_1求得加上補綴值得明文m_1
(%i20)　m1:t*(c2+2*c1-t^3)/(c2-c1+2*t^3);
(m1)　\displaystyle \frac{t(5645030747741020761470543223-t^3)}{2t^3+1632373556178779336223}

將t值帶入m_1公式
(%i21)　ev(m1,t=rhs(X[1]));
(%o21)　1234567890

將後5位補綴值刪除，得到明文M
(%i22)　M:floor(%/10^5);
(M)　12345

方法	範例
問題敘述
相同的明文m，為明文加上隨機補綴值後加密成k個密文c_i。 其中\displaystyle \alpha<\frac{k-2}{6k-3}<\frac{1}{6}，\displaystyle |\; t_i|\;\le \frac{1}{2}N^{\alpha} A_0\equiv m^3\pmod N，第0個直接三次方不加補綴 A_i\equiv (m+t_i)^3 \pmod N，第i個加上補綴後三次方 c_i\equiv A_i-A_0\equiv 3m^2t_i+3mt_i^2+t_i^3 \pmod N，第i個減第0個得到密文c_i 已知c_i,N，則透過Coppersmith方法可以回復補綴值t_i和明文m。	明文m=25152118001609140014150009191379 公鑰N=54957464841358314276864542898551 出自Cryptanalytic Attacks on RSA k=4，\displaystyle \alpha=\frac{1}{11}，\displaystyle \alpha<\frac{4-2}{6\cdot 4-3}<\frac{1}{6} t_1=761,t_2=683,t_3=714,t_4=756 其中t_1,t_2為質數，可以解出明文m 但t_3=2\cdot 3 \cdot 7\cdot 17,t_4=2^2\cdot 3^3 \cdot 7可以分解成小質數乘積，反而無法解出明文m A_0\equiv m^3 \equiv 37393323096087665763922106857101\pmod{N} A_1\equiv (m+761)^3\equiv  21904263018805857488488858542023\pmod{N} A_2\equiv (m+683)^3\equiv  37153632560891277497054197384710\pmod{N} A_3\equiv (m+714)^3\equiv  18092792378564256587234068711460\pmod{N} A_4\equiv (m+756)^3\equiv  39662861356261303674301781451309\pmod{N} 密文 c_1\equiv A_1-A_0\equiv 39468404764076506001431294583473\pmod{N} c_2\equiv A_2-A_0\equiv 54717774306161926009996633426160\pmod{N} c_3\equiv A_3-A_0\equiv 35656934123834905100176504752910\pmod{N} c_4\equiv A_4-A_0\equiv 2269538260173637910379674594208\pmod{N}
步驟1：利用恆等式建立矩陣M的右上角區塊
設下標i<j<l，定義d_{ij}=t_i t_j(t_i-t_j)，e_{ijl}=-t_it_jt_l(t_i-t_j)(t_j-t_l)(t_l-t_i) \displaystyle C_2^k個線性獨立數量d_{ij}滿足|\;d_{ij}|\;<N^{3\alpha} \displaystyle C_3^k個線性獨立數量e_{ijl}滿足|\;e_{ijl}|\;<N^{6\alpha} 滿足恆等式d_{ij}c_l+d_{jl}c_i-d_{il}c_j\equiv e_{ijl}\pmod{N} 右上角C_2^k\times C_3^k區塊為 以數對(i,j),i<j決定哪一列，以數對(i,j,l),i<j<l決定哪一行 第(i,j,l)行，第\displaystyle(i,j)=\frac{1}{2}(2k-i)(i-1)+j-i列放c_l 第(i,j,l)行，第\displaystyle(j,l)=\frac{1}{2}(2k-j)(j-1)+l-j列放c_i 第(i,j,l)行，第\displaystyle(i,l)=\frac{1}{2}(2k-i)(i-1)+l-i列放-c_j －－－－－－－－－－－－－－－－－－－－－－ 數對(i,j)在第\displaystyle \frac{1}{2}(2k-i)(i-1)+j-i列 [證明] \matrix{d_{12}&d_{13}&d_{14}&&\ldots&&d_{1k}&k-1個\cr &d_{23}&d_{24}&&&&d_{2k}&k-2個\cr &&&\ldots&&&&\cr &&&d_{i-1,i}&\ldots&&d_{i-1,k}&k-i+1個\cr &&&&d_{i,i+1}&\ldots&d_{i,j}&j-i個} (i,j)=(k-1)+(k-2)+\ldots+(k-i+1)+(j-1) \displaystyle =\frac{1}{2}(k-1+k-i+1)(i-1)+(j-i) \displaystyle =\frac{1}{2}(2k-i)(i-1)+j-i	(1,2,3)，d_{12}c_3+d_{23}c_1-d_{13}c_2\equiv e_{123}\pmod{N} (1,2,4)，d_{12}c_4+d_{24}c_1-d_{14}c_2\equiv e_{124}\pmod{N} (1,3,4)，d_{13}c_4+d_{34}c_1-d_{14}c_3\equiv e_{134}\pmod{N} (2,3,4)，d_{23}c_4+d_{34}c_2-d_{24}c_3\equiv e_{234}\pmod{N} 寫成矩陣形式 \left[\matrix{d_{12}& d_{13}& d_{14}& d_{23}& d_{24}& d_{34}}\right] \left[\matrix{c_3&c_4&&\cr -c_2&&c_4&\cr &-c_2&-c_3&\cr c_1&&&c_4\cr &c_1&&-c_3\cr &&c_1&c_2}\right]\equiv \left[\matrix{e_{123}\cr e_{124}\cr e_{134}\cr e_{234}}\right]\pmod{N} 第1行 第\displaystyle(1,2)=\frac{1}{2}(2\times4-1)(1-1)+2-1=1列放c_3 第\displaystyle(2,3)=\frac{1}{2}(2\times4-2)(2-1)+3-2=4列放c_1 第\displaystyle(1,3)=\frac{1}{2}(2\times4-1)(1-1)+3-1=2列放-c_2 第2行 第\displaystyle(1,2)=\frac{1}{2}(2\times4-1)(1-1)+2-1=1列放c_4 第\displaystyle(2,4)=\frac{1}{2}(2\times4-2)(2-1)+4-2=5列放c_1 第\displaystyle(1,4)=\frac{1}{2}(2\times4-1)(1-1)+4-1=3列放-c_2 第3行 第\displaystyle(1,3)=\frac{1}{2}(2\times4-1)(1-1)+3-1=2列放c_4 第\displaystyle(3,4)=\frac{1}{2}(2\times4-3)(3-1)+4-3=6列放c_1 第\displaystyle(1,4)=\frac{1}{2}(2\times4-1)(1-1)+4-1=3列放-c_3 第4行 第\displaystyle(2,3)=\frac{1}{2}(2\times4-2)(2-1)+3-2=4列放c_4 第\displaystyle(3,4)=\frac{1}{2}(2\times4-3)(3-1)+4-3=6列放c_2 第\displaystyle(2,4)=\frac{1}{2}(2\times4-2)(2-1)+4-2=5列放-c_3
步驟2：建立矩陣M，利用LLL找到較短的s向量
矩陣M為C_2^k+C_2^k大小的方陣， 左上角C_2^k\times C_2^k區塊為單位矩陣乘上近似N^{3\alpha}的整數值。 左下角C_3^k\times C_2^k區塊為0矩陣。 右下角C_3^k\times C_3^k區塊為單位矩陣乘上N。 右上角C_2^k\times C_3^k區塊為步驟1的係數矩陣 －－－－－－－－－－－－－－－－－－－－－－ 設列向量r前C_2^k個元是d_{ij}，後C_3^k個元是\displaystyle \frac{e_{ijl}-(d_{ij}c_l+d_{jl}c_i-d_{il}c_j)}{N} 計算rM=s 列向量s前C_2^k個元是d_{ij}N^{3\alpha}，後C_3^k個元是e_{ijl}，每個元都小於N^{6\alpha} 希望透過lattice化簡演算法找到這一列向量。 因為矩陣M為上三角矩陣，行列式值為對角線C_2^k個N^{3\alpha}和C_3^k個N相乘 \displaystyle det(M)=N^{3\alpha C_2^k+C_3^k}，因為\displaystyle \alpha<\frac{k-2}{6k-2}，所以還會大於\displaystyle  (N^{6\alpha})^{C_2^k+C_3^k} 所以s是由矩陣M個列向量所產生的較短向量，找到s的難度取決於它在短元素中的rank，若|\;t_i|\;足夠小於N^{\alpha}，然後我們希望s是整個lattice中最短元素，透過lattice化簡演算法能有效地回復s。我們在這裡不提供效率估計或成功機率，將這個方法視為啟發式攻擊(不一定100%成功)。 －－－－－－－－－－－－－－－－－－－－－－ 當\displaystyle \alpha<\frac{k-2}{6k-2}時，3\alpha C_2^k+C_3^k>6\alpha(C_2^k+C_3^k) [證明] (3\alpha C_2^k+C_3^k)-6\alpha(C_2^k+C_3^k) =(1-6\alpha)C_3^k-3\alpha C_2^k \displaystyle =(1-6\alpha)\left[\frac{k(k-1)(k-2)}{6}-\frac{k(k-1)}{2}\right] \displaystyle =\frac{k(k-1)}{2}\left[(1-6\alpha)\frac{k-2}{3}-3\alpha \right] \displaystyle \frac{k(k-1)}{2}>0為正，還需要\displaystyle (1-6\alpha)\frac{k-2}{3}-3\alpha>0亦為正 \displaystyle \frac{1-6\alpha}{\alpha}>\frac{9}{k-2} \displaystyle \frac{1}{\alpha}-6>\frac{9}{k-2} \displaystyle \frac{1}{\alpha}>\frac{9}{k-2}+\frac{6k-12}{k-2}=\frac{6k-3}{k-2} \displaystyle \alpha<\frac{k-2}{6k-3}就是當初的條件 步驟可以逆推回去	M=\left[\matrix{N^{3\alpha}&&&&&&c_3&c_4&&\cr &N^{3\alpha}&&&&&-c_2&&c_4&\cr &&N^{3\alpha}&&&&&-c_2&-c_3&\cr &&&N^{3\alpha}&&&c_1&&&c_4\cr &&&&N^{3\alpha}&&&c_1&&-c_3\cr &&&&&N^{3\alpha}&&&c_1&c_2\cr &&&&&&N&&&\cr &&0&&&&&N&&\cr &&&&&&&&N&\cr &&&&&&&&&N}\right] r=[\displaystyle d_{12},d_{13},d_{14},d_{23},d_{24},d_{34}, 　　\displaystyle \frac{e_{123}-(d_{12}c_3+d_{23}c_1-d_{13}c_2)}{N}, 　　\displaystyle \frac{e_{124}-(d_{12}c_4+d_{24}c_1-d_{14}c_2)}{N}, 　　\displaystyle \frac{e_{134}-(d_{13}c_4+d_{34}c_1-d_{14}c_3)}{N}, 　　\displaystyle \frac{e_{234}-(d_{23}c_4+d_{34}c_2-d_{24}c_3)}{N}] 計算rM=s s=[d_{12}N^{3\alpha},d_{13}N^{3\alpha},d_{14}N^{3\alpha},d_{23}N^{3\alpha},d_{24}N^{3\alpha},d_{34}N^{3\alpha}, 　　e_{123},e_{124},e_{134},e_{234}] 利用lattice化簡演算法找到短向量s，若真的找到s，從r=sM^{-1}得到r向量。
步驟3：再用Franklin和Reiter的技巧回復明文m
利用r向量計算最大公因數來回復t_i gcd(d_{1,2},d_{1,3},\ldots,d_{1,k})=gcd\{\;t_1t_2(t_1-t_2),t_1t_3(t_1-t_3),\ldots,t_1t_k(t_1-t_k) \}\; =t_1\times gcd\{\; t_2(t_1-t_2),t_3(t_1-t_3),\ldots,t_k(t_1-t_k) \}\; 希望後面的gcd足夠小而可以暴力搜尋找到t_i，再用Franklin和Reiter的技巧回復明文m。 c_1\equiv 3m^2t_1+3mt_1^2+t_1^3\pmod{N} c_2\equiv 3m^2t_2+3mt_2^2+t_2^3\pmod{N} \displaystyle m \equiv -\frac{t_1t_2^3+(c_1-t_1^3)t_2-c_2t_1}{3t_1t_2^2-3t_1^2t_2}\pmod{N}	計算gcd(r_1,r_2,r_3)得到t_1 計算gcd(r_1,r_4,r_5)得到t_2 計算gcd(r_2,r_4,r_6)得到t_3 計算gcd(r_3,r_5,r_6)得到t_4 再用Franklin和Reiter的技巧回復明文m

參考資料：
D. Coppersmith. Finding a small root of a univariate modular equation. In Proceedings of Eurocrypt 1996, volume 1070 of Lecture Notes in Computer Science, pages 155–165. Springer, 1996.
https://link.springer.com/chapter/10.1007/3-540-68339-9_14


請下載LLL.zip，解壓縮後將LLL.mac放到C:\maxima-5.46.0\share\maxima\5.46.0\share目錄下
要先載入LLL.mac才能使用LLL指令
(%i1)　load("LLL.mac");
(%o1)　C:/maxima-5.46.0/share/maxima/5.46.0/LLL.mac

公鑰N
(%i2)　N:54957464841358314276864542898551;
(N)　54957464841358314276864542898551

密文c_1,c_2,c_3,c_4
(%i3)
c:[39468404764076506001431294583473,
    54717774306161926009996633426160,
    35656934123834905100176504752910,
    2269538260173637910379674594208];
(c)　[39468404764076506001431294583473,
54717774306161926009996633426160,
35656934123834905100176504752910,
2269538260173637910379674594208]

有k=4個密文
(%i4)　k:length(c);
(k)　4

取\displaystyle \alpha=\frac{1}{11}
(%i5)　alpha:1/ceiling((6*k-3)/(k-2));
(alpha)　\displaystyle \frac{1}{11}

檢查是否符合\displaystyle \alpha<\frac{k-2}{6k-3}
(%i6)　is(alpha<(k-2)/(6*k-3));
(%o6)　true

從集合S中任選n個的組合
(%i7)
Combination(L,n):=block
([c:[],s],
if length(L)>n then
  (for r in L do
     (s:delete(r,L),
      c:unique(append(c,Combination(s,n)))
     )
  )
else
  (c:[L]),
return(c)
)$

從1\sim k中任選3個的組合
(%i8)　S:Combination(create_list(i,i,1,k),3);
(S)　[[1,2,3],[1,2,4],[1,3,4],[2,3,4]]

恆等式的係數矩陣(初始值設為0)
(%i9)　CoefficientMatrix:zeromatrix(k*(k-1)/2,k*(k-1)*(k-2)/6);
(CoefficientMatrix)　\left[\matrix{0&0&0&0\cr 0&0&0&0\cr 0&0&0&0\cr 0&0&0&0\cr 0&0&0&0\cr 0&0&0&0}\right]

恆等式的係數填入矩陣
(%i10)
for column:1 thru length(S) do
  ([i,j,l]:S[column],
   print("d",i,j,"c",l,"+d",j,l,"c",i,"-d",i,l,"c",j,"=e",i,j,l,"mod(N)"),
   dij_position: (2*k-i)*(i-1)/2+j-i,
   CoefficientMatrix[dij_position][column]:'c[l],
   djl_position: (2*k-j)*(j-1)/2+l-j,
   CoefficientMatrix[djl_position][column]:'c[ i ],
   dil_position: (2*k-i)*(i-1)/2+l-i,
   CoefficientMatrix[dil_position][column]:-'c[j]
  );
d12c3+d23c1-d13c2=e123\pmod{N}
d12c4+d24c1-d14c2=e124\pmod{N}
d13c4+d34c1-d14c3=e134\pmod{N}
d23c4+d34c2-d24c3=e234\pmod{N}
(%o10)　done

得到恆等式的係數矩陣
(%i11)　CoefficientMatrix;
(%o11)　\left[\matrix{c_3&c_4&0&0\cr -c_2&0&c_4&0\cr 0&-c_2&-c_3&0\cr c_1&0&0&c_4\cr 0&c_1&0&-c_3 \cr 0&0&c_1&c_2}\right]

單位矩陣乘上N^{3\alpha}
(%i12)　NalphaMatrix:round('N^(3*'alpha))*ident(k*(k-1)/2);
(NalphaMatrix)　\left[\matrix{ round(N^{3\alpha})&0&0&0&0&0\cr 0&round(N^{3\alpha})&0&0&0&0\cr 0&0&round(N^{3\alpha})&0&0&0\cr 0&0&0&round(N^{3\alpha})&0&0\cr 0&0&0&0&round(N^{3\alpha})&0\cr 0&0&0&0&0&round(N^{3\alpha})}\right]

單位矩陣乘上N
(%i13)　NMatrix:'N*ident(k*(k-1)*(k-2)/6);
(NMatrix)　\left[\matrix{ N&0&0&0\cr 0&N&0&0\cr 0&0&N&0\cr 0&0&0&N}\right]

零矩陣
(%i14)　ZeroMatrix:zeromatrix(k*(k-1)*(k-2)/6,k*(k-1)/2);
(ZeroMatrix)　\left[\matrix{ 0&0&0&0&0&0\cr 0&0&0&0&0&0\cr 0&0&0&0&0&0\cr 0&0&0&0&0&0}\right]

將4個子矩陣合併為矩陣M
(%i15)　M:addrow(addcol(NalphaMatrix,CoefficientMatrix),
                 addcol(ZeroMatrix,NMatrix));
(M)　\left[\matrix{ round(N^{3\alpha})&0&0&0&0&0&c_3&c_4&0&0\cr 0&round(N^{3\alpha})&0&0&0&0&-c_2&0&c_4&0\cr 0&0&round(N^{3\alpha})&0&0&0&0&-c_2&-c_3&0\cr 0&0&0&round(N^{3\alpha})&0&0&c_1&0&0&c_4\cr 0&0&0&0&round(N^{3\alpha})&0&0&c_1&0&-c_3\cr 0&0&0&0&0&round(N^{3\alpha})&0&0&c_1&c_2\cr 0&0&0&0&0&0&N&0&0&0\cr 0&0&0&0&0&0&0&N&0&0\cr 0&0&0&0&0&0&0&0&N&0\cr 0&0&0&0&0&0&0&0&0&N}\right]

將秘文c_i、公鑰N和\alpha代入矩陣M
(%i16)　M:ev(M,append(create_list('c[i ]=c[ i ],i,1,k),['N=N,'alpha=alpha]));
(M)　\left[\matrix{ 453284549&0&0&0&0&0&35656934123834905100176504752910&2269538260173637910379674594208&0&0\cr 0&453284549&0&0&0&0&-54717774306161926009996633426160&0&2269538260173637910379674594208&0\cr 0&0&453284549&0&0&0&0&-54717774306161926009996633426160&-35656934123834905100176504752910&0\cr 0&0&0&453284549&0&0&39468404764076506001431294583473&0&0&2269538260173637910379674594208\cr 0&0&0&0&453284549&0&0&39468404764076506001431294583473&0&-35656934123834905100176504752910\cr 0&0&0&0&0&453284549&0&0&39468404764076506001431294583473&54717774306161926009996633426160\cr 0&0&0&0&0&0&54957464841358314276864542898551&0&0&0\cr 0&0&0&0&0&0&0&54957464841358314276864542898551&0&0\cr 0&0&0&0&0&0&0&0&54957464841358314276864542898551&0\cr 0&0&0&0&0&0&0&0&0&54957464841358314276864542898551}\right]

用LLL進行化簡
(%i17)　B: LLL(M);
(B)　\left[\matrix{ 3062806981544531&1929302787225877&217318211327070&-1142089856961263&-2847639605402466&-1712730228981912&-7029209321862&-1864504228860&-675725901480&-5840140628952\cr 30350192335087291&-33755294424101878&-35992078252296415&10262189034662282&-12649485090783171&26238535675026248&45236596803297303&32667816209644487&17312818708413949&29899743372135217\cr -105442779364366051&-35112147780672047&15105397548793484&-162848150710657127&-47818354180209784&-39252485567286516&93004125926333727&66987628433829075&35630166688401243&61279750174710207\cr 1463179406660001&-38794361299156843&127104398505831225&15663924344526657&-53193847940963451&49664226773200368&34028664295792700&-108825608856904341&-70649517214082037&72089239243509663\cr 3555015074767161&-238087723867355568&-173132241548157991&-194306209012687083&-142531522469114792&82794656207822923&-190601615804163537&-137506222531495333&-73344681040746143&-126138894058078265\cr -195888746853505215&-146909831017677869&61301299463938392&224903490834995215&-303473180577620820&-157213305662347620&-69229181322847560&222932053809272745&145526846414244939&-148701738837887865\cr 307004212366900391&12863118098726871&599044032676087473&-259437126504168679&178098884431111932&513690572043907629&-143869831596583000&457614818638557519&299900502491764839&-303250617849528111\cr 1991487835630480536&-3276384603102473239&1317877116205579554&-212626855403276543&2241302488597860594&-3546671373603670004&4710859391231467810&-3897062285598526116&3293990032544411069&-4885714390931151907\cr 3040490414277923020&-5012515185936123906&2047825263139913024&-333778398019756331&3453617971057334454&-5468786420936785769&3167236119733748050&8002984534870019480&-15326812883649580004&-4476135453415756279\cr -10796555276199580276&17895467629929799515&-7198812808408856722&1270730505249622124&-12198909881552966272&19372609181153480002&26360591268575871440&7060175662321356&-17588138285651759645&-29785375118608432724}\right]

取第一列較短的向量s
(%i18)　s:B[1];
(s)　\left[\matrix{3062806981544531&1929302787225877&217318211327070&-1142089856961263&-2847639605402466&-1712730228981912&-7029209321862&-1864504228860&-675725901480&-5840140628952}\right]

計算r=sM^{-1}得到向量r
(%i19)　r:s.invert(M);
(r)　\left[\matrix{6756919&4256273&479430&-2519587&-6282234&-3778488&1663229&4709970&2848860&-209916}\right]

從向量r計算最大公因數得到補綴值t_i
其中t_1=761,t_2=683為質數，可以解出明文m
但t_3=714=2\cdot 3 \cdot 7\cdot 17,t_4=756=2^2\cdot 3^3 \cdot 7可以分解成小質數乘積，反而無法解出明文m
其中119僅是714的因數，126僅是756的因數
(%i21)
t:create_list(0,i,1,k)$
for i:1 thru k do
  (dij:[],
   for j:1 thru k do
      (if i<j then
         (print("取d",i,j,"=r[",index: (2*k-i-2)*(i-1)/2+j-1,"]=",r[1,index]),
          dij:append(dij,[r[1,index]])
         ),
       if i>j then
         (print("取d",j,i,"=r[",index: (2*k-j-2)*(j-1)/2+i-1,"]=",r[1,index]),
          dij:append(dij,[r[1,index]])
         )
      ),
   print("t",i,"=gcd(",dij,")=",t[ i ]:lreduce('gcd,dij))
  );
取d_{12}=r[1]=6756919
取d_{13}=r[2]=4256273
取d_{14}=r[3]=479430
t_1=gcd([6756919,4256273,479430])=761
取d_{12}=r[1]=6756919
取d_{23}=r[4]=-2519587
取d_{24}=r[5]=-6282234
t_2=gcd([6756919,-2519587,-6282234])=683
取d_{13}=r[2]=4256273
取d_{23}=r[4]=-2519587
取d_{34}=r[6]=-3778488
t_3=gcd([4256273,-2519587,-3778488])=119
取d_{14}=r[3]=479430
取d_{24}=r[5]=-6282234
取d_{34}=r[6]=-3778488
t_4=gcd([479430,-6282234,-3778488])=126
(%o21)　done　

使用Franklin和Reiter的兩個訊息有仿射關係，利用輾轉相除法找出關係式
https://math.pro/db/viewthread.php?tid=3498&page=2#pid23631
(%i22)
GCD(fx1,fx2,var):=block([temp],
fx1:expand(fx1),
fx2:expand(fx2),
while hipow(fx2,var)#1 do
  (temp:fx2,
   print(fx1,"除以",fx2,"餘式",fx2:remainder(fx1,fx2,var)),
   fx1:temp
  ),
fx2
)$

兩個多項式有相同的m
(%i24)
fx1: (m+t1)^3-m^3-c1;
fx2: (m+t2)^3-m^3-c2;
(fx1)　(t1+m)^3-m^3-c1
(fx2)　(t2+m)^3-m^3-c2

使用Franklin和Reiter輾轉相除法得到因式
(%i25)　GCD:GCD(fx1,fx2,m);
t1^3+3mt1^2+3m^2t1-c1除以t2^3+3mt2^2+3m^2t2-c2餘式\displaystyle -\frac{t1t2^3+m(3t1t2^2-3t1^2t2)+(c1-t1^3)t2-c2t1}{t2}
(GCD)　\displaystyle -\frac{t1t2^3+m(3t1t2^2-3t1^2t2)+(c1-t1^3)t2-c2t1}{t2}

得到明文m和祕文c_1,c_2和補綴值t_1,t_2的關係式
(%i26)　m:solve(GCD,m)[1];
(m)　\displaystyle m=-\frac{t1t2^3+(c1-t1^3)t2-c2t1}{3t1t2^2-3t1^2t2}

將祕文c_1,c_2和補綴值t_1,t_2的值代入關係式
(%i27)　m:ev(m,['t1=t[1],'t2=t[2],'c1=c[1],'c2=c[2]]);
(m)　\displaystyle m=-\frac{t1t2^3+(c1-t1^3)t2-c2t1}{3t1t2^2-3t1^2t2}=-\frac{14683305793124972094629922378741917}{121624542}

同餘N後得到明文m
(%i28)　m:ratsimp(rhs(m)),modulus:N;
warning: assigning 54957464841358314276864542898551, a non-prime, to 'modulus'
(m)　25152118001609140014150009191379

最後4位數為密碼
(%i29)　PinNumber:mod(m,10000);
(PinNumber)　1379

剩下的位數為訊息
(%i30)　m:floor(m/10000);
(m)　2515211800160914001415000919

以每兩位數字轉換成對應的英文字母
空白=00,A=01,B=02,\ldots,Z=26
(%i33)
Message: PinNumber$
while m#0 do
  (char:ascii(mod(m,100)+cint("A")-1),
   Message:concat(char,Message),
   m:floor(m/100)
  )$
Message;
(%o33)　YOUR@PIN@NO@IS1379

將@取代為空白
(%i34)　Message:ssubst(" ","@",Message);
(Message)　YOUR PIN NO IS1379

3-1.針對大的r，分解RSA公鑰n=p^rq
在RSA加密演算法中，為了加快解密速度，採用Quisquater-Couvreur方法，先各別計算m_p\equiv c^{d_p}\pmod{p}和m_q\equiv c^{d_q}\pmod{q}，再利用中國餘數定理將m_p和m_q組合出明文m，範例請見https://math.pro/db/viewthread.php?tid=1500&page=1#pid7268。
若公鑰e取較小的數，但仍要考慮到低次方攻擊，加密過程可以減少計算而且快速，另一方面，解密的私鑰d的位元數必須要比公鑰n位元數的\displaystyle \frac{1}{4}還大以避免Wiener攻擊，所以解密過程變成是整個RSA中計算量最大的。
所以Takagi提出另一種RSA類型的加密演算法，將公鑰n從原本的pq改為p^rq，以減少解密過程時的計算量，以下是原本的RSA加密演算法和Takagi所提出的RSA方案的比較表。

RSA　n=pq	RSA　n=p^rq
產生公鑰和私鑰
產生兩個隨機質數p,q，n=pq 計算\phi(n)=(p-1)(q-1) 找e,d符合ed\equiv 1 \pmod{\phi(n)}和GCD(e,p)=1 e,n是公鑰，d,p,q是私鑰	產生兩個隨機質數p,q，n=p^r q 計算L=LCM(p-1.q-1) 找e,d符合ed \equiv 1\pmod{L}和GCD(e,p)=1 e,n是公鑰，d,p,q是私鑰
加密
c\equiv m^e \pmod{n} m為明文，c為密文	相同
解密
m\equiv c^d \pmod{n}	無法使用m\equiv c^d \pmod{n}來回復明文
利用中國餘數定理快速解密
m\equiv c^d \pmod{pq} 分解成m_p\equiv c^d \pmod{p}，m_q\equiv c^d\pmod{q} 設d_p\equiv d\pmod{p-1}，m_p\equiv c^{d_p}\pmod{p} 設d_q\equiv d\pmod{q-1}，m_q\equiv c^{d_q}\pmod{q} 利用中國餘數定理 m\equiv \left(m_p\cdot q(q^{-1}\pmod{p})+m_q\cdot p(p^{-1}\pmod{q})\right)\pmod{n}	m_q\equiv c^d\pmod{q} 設d_q\equiv d\pmod{q-1}，m_q\equiv c^{d_q}\pmod{q} 利用中國餘數定理 m\equiv (m_p\cdot q(q^{-1}\pmod{p^r})+m_q\cdot p^r(p^{-r}\pmod{q}))\pmod{n} m_p\equiv m\pmod{p^r}的m_p無法用m_p\equiv c^d\pmod{p^r}直接計算，改用以下的方法計算m_p

假設密文c經p^r同餘後為c_p，明文m經p^r同餘後為m_p，密文c_p和明文m_p的關係為c_p\equiv m_p^e \pmod{p^r}。
設明文m_p以p進位表示成
m_p\equiv K_0+pK_1+p^2K_2+\ldots+p^{r-1}K_{r-1}\pmod{p^r}。
代入c_p\equiv m_p^e\pmod{p^r}
得到c_p\equiv (K_0+pK_1+p^2K_2+\ldots+p^{r-1}K_{r-1})^e\pmod{p^r}
設關係式
c_p\equiv (\bbox[border:1px solid black]{K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1}}+\bbox[border:1px solid black]{p^iK_i})^e\pmod{p^{i+1}}，i=1,2,\ldots,r-1
利用二項式定理展開
c_p\equiv C_0^e(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^e(p^iK_i)^0
　　+C_1^e(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^{e-1}(p^iK_i)^1
　　+C_2^e(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^{e-2}(p^iK_i)^2　超過p^{i+1}
　　　　　...
　　+C_e^e(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^{0}(p^iK_i)^e超過p^{i+1}\pmod{p^{i+1}}
化簡後
c_p\equiv (K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^e+e(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^{e-1}p^iK_i\pmod{p^{i+1}}
定義A_i=K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1}，
　　F_i=(A_i)^e\pmod{p^{i+1}}，
　　(K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1})^{e-1}=F_iA_i^{-1}
關係式為c_p\equiv F_i+eF_iA_i^{-1}p^iK_{i}\pmod{p^{i+1}}，i=1,2,\ldots,r-1
移項可得\displaystyle K_{i}\equiv \frac{(c_p-F_i)(eF_i)^{-1}A_i}{p^i}\pmod{p^{i+1}}，i=1,2,\ldots,r-1
就可以計算出K_1,K_2,\ldots,K_{r-1}的值，K_0\equiv c^{d_p}\pmod{p}，進而求出m_p=K_0+pK_1+p^2K_2+\ldots+p^{r-1}K_{r-1}。

範例的p=61,q=53,e=17出自
https://en.wikipedia.org/wiki/RSA_(cryptosystem)#Example

設r=3
公鑰n=p^rq=61^3\cdot 53=12029993
L=LCM(p-1,q-1)=LCM(60,52)=780
私鑰d\equiv e^{-1}\pmod{L}\equiv 17^{-1}\equiv 413\pmod{780}
設明文m=1234567
密文c\equiv m^e\pmod{n}
　　　\equiv 1234567^{17}\pmod{12029993}
　　　\equiv 10259995

其他私鑰d方案

L=LCM(p-1,q-1) d\equiv e^{-1}\pmod{L}	\phi(n)=p^{r-1}(p-1)(q-1) d\equiv e^{-1}\pmod{\phi(n)}	\lambda(n)=LCM(p^{r-1}(p-1),(q-1)) d\equiv e^{-1}\pmod{\lambda(n)}
L=LCM(60,52)=780 d\equiv 17^{-1}\equiv 413\pmod{L} 私鑰d較小，但不能直接計算m\equiv c^d \pmod{n}來回復明文，改用以下的中國餘數定理來加速解密。	\phi(n)=61^2\cdot 60\cdot 52=11609520 d\equiv 17^{-1}\equiv 682913\pmod{\phi(n)} 雖然可以直接計算m\equiv c^d \pmod{n}來回復明文， 但解密速度會比n=pq且使用中國餘數定理的RSA還慢，變成使用n=p^rq的RSA沒有顯著的優點。	\lambda(n)=LCM(61^2\cdot 60,52)=2902380 d\equiv 17^{-1}\equiv 682913\pmod{\lambda(n)} 雖然可以直接計算m\equiv c^d \pmod{n}來回復明文， 但解密速度會比n=pq且使用中國餘數定理的RSA還慢，變成使用n=p^rq的RSA沒有顯著的優點。

使用中國餘數定理解密。

虛擬碼

範例

A_i=K_0+pK_1+p^2K_2+\ldots+p^{i-1}K_{i-1}改成遞迴關係式A_{i+1}=A_i+p^iK_i F_i=(A_i)^e\pmod{p^{i+1}} 計算\displaystyle K_{i}\equiv \frac{(c_p-F_i)(eF_i)^{-1}A_i}{p^i}\pmod{p^{i+1}}，i=1,2,\ldots,r-1 分解成 E_i=c-F_i\pmod{p^{i+1}} \displaystyle B_i=\frac{E_i}{p^i} K_i=[(eF_i)^{-1}A_iB_i]\pmod{p} －－－－－－－－－－－－－－－－ 解密 輸入：d,p,q,e,r,c 輸出：m (1)　d_p=d\pmod{p-1} 　　d_q=d\pmod{q-1} (2)　K_0=c^{d_p}\pmod{p} 　　m_q=c^{d_q}\pmod{q} (3)　A_1=K_0 　　for i=1 to (r-1) do 　　　F_i=A_i^e \pmod{p^{i+1}} 　　　E_i=(c-F_i)\pmod{p^{i+1}} 　　　\displaystyle B_i=\frac{E_i}{p^i}為整數 　　　K_i=[(eF_i)^{-1}A_iB_i]\pmod{p} 　　　A_{i+1}=A_i+p^iK_i為整數 (4)　m_p=A_r (5)　p_1=[(p^r)^{-1}]\pmod{q} 　　q_1=[q^{-1}]\pmod{p^r} (6)　m=[q_1qm_p+p_1p^rm_q]\pmod{p^rq}

輸入：d=413，p=61，q=53，e=17，r=3，c=10259995 輸出：m=1234567 (1)d_p\equiv 413\pmod{60}\equiv 53 　d_q\equiv 413\pmod{52}\equiv 49 (2)K_0\equiv 10259995^{53}\equiv 49\pmod{61} 　m_q\equiv 10259995^{49}\equiv 38\pmod{53} (3)A_1=49 當i=1時 　F_1=A_1^e= 49^{17}  = 527 \pmod{61 ^2} 　E_1=(c-F_1)=(10259995-527)=671 \pmod{61 ^2} 　\displaystyle B_1=\frac{E_1}{p}=\frac{671}{61}=11 　K_1=(eF_1)^{-1} A_1B_1=( 17 \cdot 527 )^{-1} 49 \cdot 11 = 47 \pmod{61} 　A_2=A_1+p \cdot K_1= 49 + 61  \cdot 47 = 2916 當i=2時 　F_2=A_2^e= 2916^{17}= 57013 \pmod{61^3} 　E_2=(c-F_2)=(10259995-57013)= 215818 \pmod{61 ^3} 　\displaystyle B_2=\frac{E_2}{p^2}=\frac{215818}{3721}= 58 　K_2=(eF_2)^{-1} A_2B_2=( 17 \cdot 57013 )^{-1} 2916 \cdot 58 = 26 \pmod{61} 　A_3=A_2+p^2 \cdot K_2= 2916 + 61 ^2 \cdot 26 = 99662 (4)m_p=A_r=99662 (5)p_1\equiv 61^{-3}\equiv 50\pmod{53} 　q_1\equiv 53^{-1}\equiv 12848 \pmod{61^3} (6)m\equiv [12848\cdot 53\cdot 99662+50\cdot 61^3 \cdot 38]\pmod{61^3\cdot 53} 　　\equiv 1234567 \pmod{61^3\cdot 53}

註：
1.論文的A[0]index從0開始，但maxima的list從1開始，所以將A_{i-1}改為A_i和A_i改為A_{i+1}和A_{r-1}改為A_r。
2.原論文公鑰n=p^kq，順應分解RSA公鑰n=p^rq論文，符號k另有他用，p的次方改為r。

參考資料：
Takagi, T. (1998). Fast RSA-type cryptosystem modulo p k q . In: Krawczyk, H. (eds) Advances in Cryptology — CRYPTO '98. CRYPTO 1998. Lecture Notes in Computer Science, vol 1462. Springer, Berlin, Heidelberg. https://doi.org/10.1007/BFb0055738
https://link.springer.com/chapter/10.1007/bfb0055738
https://link.springer.com/conten ... f?pdf=inline%20link



選兩個不同的質數
(%i2)
p:61;
q:53;
(p)　61
(q)　53

公鑰n=p^rq
(%i4)
r:3;
n:p^r*q;
(r)　3
(n)　12029993

公鑰e
(%i5)　e:17;
(e)　17

私鑰d
(%i6)　d:inv_mod(e,lcm(p-1,q-1));
(d)　413

明文m
(%i7)　m:1234567;
(m)　1234567

用m^e\pmod{n}進行加密，密文為c
(%i8)　c:power_mod(m,e,n);
(c)　10259995

用c^d\pmod{n}直接計算m\equiv c^d\pmod(n)無法得到原訊息m
(%i9)　m:power_mod(c,d,n);
(m)　2330554

利用中國餘數定理才能解密
各變數初始化
(%i18)
dp:mod(d,p-1);
dq:mod(d,q-1);
K:create_list(0,i,1,r);
K[1]:power_mod(c,dp,p);
mq:power_mod(c,dq,q);
A:create_list(0,i,1,r);
B:create_list(0,i,1,r);
F:create_list(0,i,1,r);
A[1]:K[1];
(dp)　53
(dq)　49
(K)　[0,0,0]
(K[1])　49
(mq)　38
(A)　[0,0,0]
(B)　[0,0,0]
(F)　[0,0,0]
(A[1])　49

利用F_i=(A_i)^e，\displaystyle K_i\equiv \frac{(c-F_i)(eF_i)^{-1}A_i}{p^i}\pmod{p^{i+1}}，i=1,2,\ldots,r-1，求K_i值
(%i9)
for i:1 thru r-1 do
  (print("F[",i,"]=A[",i,"]"^"e=",A[ i ],""^e," =",F[ i ]:power_mod(A[ i ],e,p^(i+1)),"(mod",p,""^(i+1),")"),
   print("E[",i,"]=(c-F[",i,"])=(",c,"-",F[ i ],")=",E[ i ]:mod(c-F[ i ],p^(i+1)),"(mod",p,""^(i+1),")"),
    print("B[",i,"]=E[",i,"]/p"^i,"=",E[ i ],"/",p^i,"=",B[ i ]:E[ i ]/p^i),
    print("K[",i,"]=(eF[",i,"])"^"-1","A[",i,"]B[",i,"]=(",e,"*",F[ i ],")"^"-1",A[ i ],"*",B[ i ],"=",
            K[ i ]:mod(inv_mod(e*F[ i ],p)*A[ i ]*B[ i ],p),"(mod",p,")"),
    print("A[",i+1,"]=A[",i,"]+p"^i,"*K[",i,"]=",A[ i ],"+",p,""^i,"*",K[ i ],"=",A[i+1]:A[ i ]+p^i*K[ i ]),
    print("－－－－－－")
  )$
F[1]=A[1]^e=49^{17}=527\pmod{61^2}
E[1]=(c-F[1])=(10259995-527)=671\pmod{61^2}
B[1]=E[1]/p=671/61=11
K[1]=(eF[1])^{-1}A[1]B[1]=(17*527)^{-1}49*11=47\pmod{61}
A[2]=A[1]+p*K[1]=49+61*47=2916
－－－－－－
F[2]=A[2]^e=2916^{17}=57013\pmod{61^3}
E[2]=(c-F[2])=(10259995-57013)=215818\pmod{61^3}
B[2]=E[2]/p^2=215818/3721=58
K[2]=(eF[2])^{-1}A[2]B[2]=(17*57013)^{-1}2916*58=26\pmod{61}
A[3]=A[2]+p^2*K[2]=2916+61^2*26=99662
－－－－－－

得到明文m
(%i23)
mp:A[r];
p1:inv_mod(p^r,q);
q1:inv_mod(q,p^r);
m:mod(q1*q*mp+p1*p^r*mq,p^r*q);
(mp)　99662
(p1)　50
(q1)　12848
(m)　1234567

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