用maxima學密碼學－Lattice Reduction

4.改進LLL方法－遞迴LLL

原本LLL方法是針對整個Lattice進行LLL運算，而遞迴LLL一開始先對前兩個向量進行LLL運算，之後再加入第三個向量進行LLL運算...，每次都新增一個向量進行LLL運算直到整個Lattice都化簡為止。

範例出自https://thomas-plantard.github.io/pdf/talk10-scn.pdf

這個LLL副程式略有修改，參數\( m \)代表將第1~m列向量進行LLL運算。
(%i1)
LLL(v,m):=block
  (vstar:copymatrix(v),
   delta:3/4,
   mu:zeromatrix(m,m),
   for i:1 thru m do
     (for j:1 thru i do
        (mu[ i ][j]:v[ i ].vstar[j]/vstar[j].vstar[j],
         vstar[ i ]:v[ i ]-sum(mu[ i ][j]*vstar[j],j,1,i-1)
        )
     ),
   k:2,
   while k<=m do
     (for i:k-1 thru 1 step -1 do
        (q:round(mu[k][ i ]),
         if q#0 then
           (v:rowop(v,k,i,q),
            mu:rowop(mu,k,i,q)
           )
        ),
      if (vstar[k].vstar[k])>=(delta-mu[k][k-1]^2)*(vstar[k-1].vstar[k-1]) then
        (k:k+1
        )
      else
        (v:rowswap(v,k,k-1),
         vstar[1]:v[1],
         for i:1 thru m do
           (for j:1 thru i do
             (mu[ i ][j]:v[ i ].vstar[j]/vstar[j].vstar[j],
              vstar[ i ]:v[ i ]-sum(mu[ i ][j]*vstar[j],j,1,i-1)
             )
           ),
         if k>2 then
           (k:k-1
           )
        )
     ),
  return(v)
)$

(%i2)
B:matrix([1031,0,0,0,0],
         [354,1,0,0,0],
         [322,0,1,0,0],
         [916,0,0,1,0],
         [426,0,0,0,1]);
(%o2)　\( \left[ \matrix{1031 & 0 & 0 & 0 & 0 \cr 354 & 1 & 0 & 0 & 0 \cr 322 & 0 & 1 & 0 & 0 \cr 916 & 0 & 0 & 1 & 0 \cr 426 & 0 & 0 & 0 & 1} \right] \)

(%i3)
for i:2 thru length(B) do
  (B: LLL(B,i),
   print(B)
  );
\( \left[ \matrix{-31 & -3 & 0 & 0 & 0 \cr 13 & -32 & 0 & 0 & 0 \cr 322 & 0 & 1 & 0 & 0 \cr 916 & 0 & 0 & 1 & 0 \cr 426 & 0 & 0 & 0 & 1} \right] \)

\( \left[ \matrix{-1 & 2 & 1 & 0 & 0 \cr 0 & -6 & 13 & 0 & 0 \cr -27 & -11 & -4 & 0 & 0 \cr 916 & 0 & 0 & 1 & 0 \cr 426 & 0 & 0 & 0 & 1} \right] \)

\( \left[ \matrix{-1 & 2 & 1 & 0 & 0 \cr -4 & -3 & 4 & 2 & 0 \cr 4 & 0 & 5 & 5 & 0 \cr 0 & -3 & 4 & -7 & 0 \cr 426 & 0 & 0 & 0 & 1} \right] \)

\( \left[ \matrix{-1 & 2 & 1 & 0 & 0 \cr 0 & -2 & 2 & 1 & -2 \cr -2 & -2 & 3 & -3 & 1 \cr -2 & 1 & -1 & 4 & 1 \cr 3 & 0 & 3 & 0 & 5} \right] \)
(%o3)　done

另一個範例出自http://www1.uwindsor.ca/sac2012/system/files/2a_Plantard.pdf
(%i4)
B:matrix([86670401,0,0,0,0,0,0,0],
         [38009011,1,0,0,0,0,0,0],
         [10117311,0,1,0,0,0,0,0],
         [38269415,0,0,1,0,0,0,0],
         [45874978,0,0,0,1,0,0,0],
         [33538152,0,0,0,0,1,0,0],
         [61611560,0,0,0,0,0,1,0],
         [66174289,0,0,0,0,0,0,1]);
(%o4)　\( \left[\matrix{
86670401&0&0&0&0&0&0&0 \cr
38009011&1&0&0&0&0&0&0 \cr
10117311&0&1&0&0&0&0&0 \cr
38269415&0&0&1&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

(%i5)
for i:2 thru length(B) do
  (B: LLL(B,i),
   print(B)
  );
\( \left[\matrix{
-3227&-3165&0&0&0&0&0&0 \cr
-14111&13018&0&0&0&0&0&0 \cr
10117311&0&1&0&0&0&0&0 \cr
38269415&0&0&1&0&0&0&0 \cr
45874979&0&0&0&1&0&0&0 \cr
33538154&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
-24&153&-215&0&0&0&0&0 \cr
183&242&76&0&0&0&0&0 \cr
-920&440&343&0&0&0&0&0 \cr
38266415&0&0&1&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
8&-27&-66&42&0&0&0&0 \cr
-40&-45&-47&-38&0&0&0&0 \cr
0&126&-18&-23&0&0&0&0 \cr
103&26&0&-53&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
-17&-31&-5&6&1&0&0&0 \cr
24&-20&-6&4&7&0&0&0 \cr
-3&-7&-45&3&17&0&0&0 \cr
8&4&-13&-14&-36&0&0&0 \cr
13&0&15&-35&24&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
-7&-6&-14&11&4&-7&0&0 \cr
-14&-9&-3&-1&-15&-6&0&0 \cr
-2&15&14&10&-6&4&0&0 \cr
8&-14&5&13&-14&-2&0&0 \cr
-5&11&-6&-10&-12&12&0&0 \cr
4&-16&12&-4&12&13&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
5&-4&-8&0&1&10&-1&0 \cr
9&-1&4&-6&-7&-8&2&0 \cr
1&-4&0&6&-12&-1&4&0 \cr
4&8&-4&9&4&0&-3&0 \cr
3&-2&-11&-4&-5&-3&-6&0 \cr
4&-9&-7&9&3&-2&7&0 \cr
7&-10&5&7&-2&-1&-4&0 \cr
66174289&0&0&0&0&0&0&1 }\right] \)

\( \left[\matrix{
-2&6&-1&-8&0&1&-3&1 \cr
2&-3&-8&1&3&-1&4&1 \cr
\color{red}{3}&\color{red}{0}&\color{red}{-3}&\color{red}{-5}&\color{red}{-3}&\color{red}{2}&\color{red}{5}&\color{red}{5} \cr
-2&-1&-6&4&-6&-3&0&-3 \cr
5&-4&4&-1&-2&0&-7&1 \cr
4&2&3&-1&-4&1&-1&-7 \cr
-2&3&-4&0&0&11&2&-2 \cr
5&11&1&3&-2&3&-2&4 }\right] \)
(%o5)　done




執行結果第三列向量和pdf檔所寫的不一樣，我還以為我的LLL程式出了問題，我還另外用maple再寫一次來測試，結果是相同的。

maple檔下載遞迴LLL1.rar
\( \left[\matrix{
-2&6&-1&-8&0&1&-3&1 \cr
2&-3&-8&1&3&-1&4&1 \cr
\color{red}{3}&\color{red}{0}&\color{red}{-3}&\color{red}{-5}&\color{red}{-3}&\color{red}{2}&\color{red}{5}&\color{red}{5} \cr
-2&-1&-6&4&-6&-3&0&-3 \cr
5&-4&4&-1&-2&0&-7&1 \cr
4&2&3&-1&-4&1&-1&-7 \cr
-2&3&-4&0&0&11&2&-2 \cr
5&11&1&3&-2&3&-2&4 }\right] \)\( \left[\matrix{
-2&6&-1&-8&0&1&-3&1 \cr
2&-3&-8&1&3&-1&4&1 \cr
\color{red}{3}&\color{red}{-3}&\color{red}{6}&\color{red}{2}&\color{red}{-6}&\color{red}{2}&\color{red}{4}&\color{red}{3} \cr
-2&-1&-6&4&-6&-3&0&-3 \cr
5&-4&4&-1&-2&0&-7&1 \cr
4&2&3&-1&-4&1&-1&-7 \cr
-2&3&-4&0&0&11&2&-2 \cr
5&11&1&3&-2&3&-2&4 }\right] \)





＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
在 http://www1.uwindsor.ca/sac2012/system/files/2a_Plantard.pdf 還提到另外一種遞迴LLL方法，先將兩個向量為一組進行LLL運算，再以四個向量為一組進行LLL運算，再以八個向量為一組進行LLL運算…，直到全部的Lattice向量都化簡完為止。

這個LLL副程式略有修改，參數\( m1,m2 \)代表將第m1~m2列向量進行LLL運算。
(%i1)
LLL(v,m1,m2):=block
  (vstar:copymatrix(v),
   delta:3/4,
   m:length(v),
   mu:zeromatrix(m,m),
   for i:m1 thru m2 do
     (for j:m1 thru i do
        (mu[ i ][j]:v[ i ].vstar[j]/vstar[j].vstar[j],
         vstar[ i ]:v[ i ]-sum(mu[ i ][j]*vstar[j],j,1,i-1)
        )
     ),
   k:m1+1,
   while k<=m2 do
     (for i:k-1 thru 1 step -1 do
        (q:round(mu[k][ i ]),
         if q#0 then
           (v:rowop(v,k,i,q),
            mu:rowop(mu,k,i,q)
           )
        ),
      if (vstar[k].vstar[k])>=(delta-mu[k][k-1]^2)*(vstar[k-1].vstar[k-1]) then
        (k:k+1
        )
      else
        (v:rowswap(v,k,k-1),
         vstar[m1]:v[m1],
         for i:m1 thru m2 do
           (for j:m1 thru i do
             (mu[ i ][j]:v[ i ].vstar[j]/vstar[j].vstar[j],
              vstar[ i ]:v[ i ]-sum(mu[ i ][j]*vstar[j],j,1,i-1)
             )
           ),
         if k>m1+1 then
           (k:k-1
           )
        )
     ),
  return(v)
)$

(%i2)
B:matrix([86670401,0,0,0,0,0,0,0],
         [38009011,1,0,0,0,0,0,0],
         [10117311,0,1,0,0,0,0,0],
         [38269415,0,0,1,0,0,0,0],
         [45874978,0,0,0,1,0,0,0],
         [33538152,0,0,0,0,1,0,0],
         [61611560,0,0,0,0,0,1,0],
         [66174289,0,0,0,0,0,0,1]);
(%o2)　\( \left[ \matrix{
86670401&0&0&0&0&0&0&0 \cr
38009011&1&0&0&0&0&0&0 \cr
10117311&0&1&0&0&0&0&0 \cr
38269415&0&0&1&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 } \right] \)

(%i3)
d:length(B)$
for i:1 thru ceiling(log(d)/log(2)) do
  (for j:1 thru ceiling(d/2^i) do
     (m1: (j-1)*2^i+1,
      m2:j*2^i,
      if m2>d then
        (m2:d),
      B: LLL(B,m1,m2),
      print("第",m1,"到",m2,"列向量進行LLL運算"),
      print(B)
     )
  );
第1到2列向量進行LLL運算
\( \left[ \matrix{
-3227&-3165&0&0&0&0&0&0 \cr
-14111&13018&0&0&0&0&0&0 \cr
10117311&0&1&0&0&0&0&0 \cr
38269415&0&0&1&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 } \right] \)

第3到4列向量進行LLL運算
\( \left[ \matrix{
-3227&-3165&0&0&0&0&0&0 \cr
-14111&13018&0&0&0&0&0&0 \cr
1391&0&4036&-1067&0&0&0&0 \cr
-8121&0&3949&-1044&0&0&0&0 \cr
45874978&0&0&0&1&0&0&0 \cr
33538152&0&0&0&0&1&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 } \right] \)

第5到6列向量進行LLL運算
\( \left[ \matrix{
-3227&-3165&0&0&0&0&0&0 \cr
-14111&13018&0&0&0&0&0&0 \cr
1391&0&4036&-1067&0&0&0&0 \cr
-8121&0&3949&-1044&0&0&0&0 \cr
1348&0&0&0&2830&-3871&0&0 \cr
10894&0&0&0&-2009&2748&0&0 \cr
61611560&0&0&0&0&0&1&0 \cr
66174289&0&0&0&0&0&0&1 } \right] \)

第7到8列向量進行LLL運算
\( \left[ \matrix{
-3227&-3165&0&0&0&0&0&0 \cr
-14111&13018&0&0&0&0&0&0 \cr
1391&0&4036&-1067&0&0&0&0 \cr
-8121&0&3949&-1044&0&0&0&0 \cr
1348&0&0&0&2830&-3871&0&0 \cr
10894&0&0&0&-2009&2748&0&0 \cr
3&0&0&0&0&0&2248&-2093 \cr
29437&0&0&0&0&0&29&-27 } \right] \)

第1到4列向量進行LLL運算
\( \left[ \matrix{
8&-27&-66&42&0&0&0&0 \cr
-40&-45&-47&-38&0&0&0&0 \cr
0&126&-18&-23&0&0&0&0 \cr
\color{red}{-103}&\color{red}{-26}&0&\color{red}{53}&0&0&0&0 \cr
1348&0&0&0&2830&-3871&0&0 \cr
10894&0&0&0&-2009&2748&0&0 \cr
3&0&0&0&0&0&2248&-2093 \cr
29437&0&0&0&0&0&29&-27 } \right] \)

第5到8列向量進行LLL運算
\( \left[ \matrix{
8&-27&-66&42&0&0&0&0 \cr
-40&-45&-47&-38&0&0&0&0 \cr
0&126&-18&-23&0&0&0&0 \cr
-103&-26&0&53&0&0&0&0 \cr
-22&0&0&0&-7&-19&-36&48 \cr
\color{red}{-93}&0&0&0&\color{red}{25}&\color{red}{2}&\color{red}{5}&\color{red}{-23} \cr
\color{red}{32}&0&0&0&\color{red}{97}&\color{red}{-13}&\color{red}{-63}&\color{red}{-2} \cr
\color{red}{-1}&0&0&0&\color{red}{25}&\color{red}{111}&\color{red}{-93}&\color{red}{13} } \right] \)

第1到8列向量進行LLL運算
\( \left[ \matrix{
\color{red}{-2}&\color{red}{3}&\color{red}{8}&\color{red}{-1}&\color{red}{-3}&\color{red}{1}&\color{red}{-4}&\color{red}{-1} \cr
\color{red}{-3}&\color{red}{0}&\color{red}{3}&\color{red}{5}&\color{red}{3}&\color{red}{-2}&\color{red}{-5}&\color{red}{-5} \cr
\color{red}{2}&\color{red}{1}&\color{red}{6}&\color{red}{-4}&\color{red}{6}&\color{red}{3}&\color{red}{0}&\color{red}{3} \cr
\color{red}{-4}&\color{red}{-2}&\color{red}{-3}&\color{red}{1}&\color{red}{4}&\color{red}{-1}&\color{red}{1}&\color{red}{7} \cr
\color{red}{-2}&\color{red}{6}&\color{red}{-1}&\color{red}{-8}&\color{red}{0}&\color{red}{1}&\color{red}{-3}&\color{red}{1} \cr
\color{red}{-5}&\color{red}{4}&\color{red}{-4}&\color{red}{1}&\color{red}{2}&\color{red}{0}&\color{red}{7}&\color{red}{-1} \cr
\color{red}{3}&\color{red}{-1}&\color{red}{0}&\color{red}{-1}&\color{red}{-2}&\color{red}{11}&\color{red}{-5}&\color{red}{-1} \cr
5&11&1&3&-2&3&-2&4 } \right] \)
(%o4)　done




紅色數字都是執行結果和pdf檔所寫的不一樣的地方，只是再用maple測試，執行結果和我的相同，或許原作者所使用的LLL程式多了其他的設定(\( \delta=0.75 \)改為\( \delta=0.99 \))

maple下載遞迴LLL2.rar
\( \left[ \matrix{
\color{red}{-2}&\color{red}{3}&\color{red}{8}&\color{red}{-1}&\color{red}{-3}&\color{red}{1}&\color{red}{-4}&\color{red}{-1} \cr
\color{red}{-3}&\color{red}{0}&\color{red}{3}&\color{red}{5}&\color{red}{3}&\color{red}{-2}&\color{red}{-5}&\color{red}{-5} \cr
\color{red}{2}&\color{red}{1}&\color{red}{6}&\color{red}{-4}&\color{red}{6}&\color{red}{3}&\color{red}{0}&\color{red}{3} \cr
\color{red}{-4}&\color{red}{-2}&\color{red}{-3}&\color{red}{1}&\color{red}{4}&\color{red}{-1}&\color{red}{1}&\color{red}{7} \cr
\color{red}{-2}&\color{red}{6}&\color{red}{-1}&\color{red}{-8}&\color{red}{0}&\color{red}{1}&\color{red}{-3}&\color{red}{1} \cr
\color{red}{-5}&\color{red}{4}&\color{red}{-4}&\color{red}{1}&\color{red}{2}&\color{red}{0}&\color{red}{7}&\color{red}{-1} \cr
\color{red}{3}&\color{red}{-1}&\color{red}{0}&\color{red}{-1}&\color{red}{-2}&\color{red}{11}&\color{red}{-5}&\color{red}{-1} \cr
5&11&1&3&-2&3&-2&4 } \right] \)\( \left[ \matrix{
\color{red}{-2}&\color{red}{6}&\color{red}{-1}&\color{red}{-8}&\color{red}{0}&\color{red}{1}&\color{red}{-3}&\color{red}{1} \cr
\color{red}{2}&\color{red}{-3}&\color{red}{-8}&\color{red}{1}&\color{red}{3}&\color{red}{-1}&\color{red}{4}&\color{red}{1} \cr
\color{red}{3}&\color{red}{-3}&\color{red}{6}&\color{red}{2}&\color{red}{-6}&\color{red}{2}&\color{red}{4}&\color{red}{3} \cr
\color{red}{-2}&\color{red}{-1}&\color{red}{-6}&\color{red}{4}&\color{red}{-6}&\color{red}{-3}&\color{red}{0}&\color{red}{-3} \cr
\color{red}{5}&\color{red}{-4}&\color{red}{4}&\color{red}{-1}&\color{red}{-2}&\color{red}{0}&\color{red}{-7}&\color{red}{1} \cr
\color{red}{4}&\color{red}{2}&\color{red}{3}&\color{red}{-1}&\color{red}{-4}&\color{red}{1}&\color{red}{-1}&\color{red}{-7} \cr
\color{red}{-2}&\color{red}{3}&\color{red}{-4}&\color{red}{0}&\color{red}{0}&\color{red}{11}&\color{red}{2}&\color{red}{-2} \cr
5&11&1&3&-2&3&-2&4 } \right] \)



參考資料：
Recursive Lattice Reduction，Thomas Plantard, Willy Susilo
https://thomas-plantard.github.io/pdf/PlaSus10.pdf
簡報檔
https://thomas-plantard.github.io/pdf/talk10-scn.pdf

Lattice Reduction for Modular Knapsack，Thomas Plantard, Willy Susilo, Zhenfei Zhang
https://thomas-plantard.github.io/pdf/PlaSusZha12.pdf
簡報檔
http://www1.uwindsor.ca/sac2012/system/files/2a_Plantard.pdf

5.改進LLL方法－PotLLL

原本的LLL只考慮相鄰的兩列向量( \( v_i^*、v_{i-1}^* \) )是否滿足Lovász condition( \( \displaystyle \Vert\; v_i^* \Vert\;^2 \ge (\frac{3}{4}-\mu_{i,i-1}^2)\Vert\; v_{i-1}^* \Vert\;^2 \) )，違反的話則將\( v_{k-1},v_k \)交換。而PotLLL考慮是否滿足(\( \delta \cdot Pot(B)\le Pot(\sigma_{k,l}B) \))，違反的話則將第\( l \)列向量插入第\( k \)列向量的位置。


\( \sigma_{k,l}(i)=\cases{i \cr l \cr i-1}\matrix{\hbox{for　i<k or i>l}, \cr \hbox{for　i=k,　} \cr \hbox{for　k<i ≦ l.}} \)代表將第\( l \)列向量插入第\(k\)列向量。
\( B \)是論文用來代表lattice，我的maxima程式用\( v \)表示。
\( \displaystyle Pot(B):= \prod_{i=1}^n vol(L(b_1,\ldots,b_i))^2= \prod_{i=1}^n \Vert\; b_i^* \Vert\;^{2(n-i+1)} \)計算lattice的Potential值
以上式子就是將lattice \(B\)分成很多sub lattice
\( B_1=L(b_1) \)，由\( b_1 \)向量所組成的lattice，其體積的平方為\( \Vert\; b_1^* \Vert\;^2 \)
\( B_2=L(b_1,b_2) \)，由\( b_1,b_2 \)向量組成的lattice，其體積的平方為\( \Vert\; b_1^* \Vert\;^2 \cdot \Vert\; b_2^* \Vert\;^2 \)
\( \ldots \)
\( B_n=L(b_1,b_2,\ldots,b_n) \)，由\( b_1,b_2,\ldots,b_n \)向量組成的lattice，其體積的平方為\( \Vert\; b_1^* \Vert\;^2 \cdot \Vert\; b_2^* \Vert\;^2 \ldots \Vert\; b_n^* \Vert\;^2 \)
將全部sub lattice體積的平方全部乘起來\( \displaystyle \prod_{i=1}^n \Vert\; b_i^* \Vert\;^{2(n-i+1)} \)就是Potential值


以\( B=\left[ \matrix{9&2&7 \cr 8&6&1 \cr 3&2&6 } \right] \)為例，計算\( B \)的Potential值
先計算GramSchmit正交化後的向量\( b^*=\left[ \matrix{\displaystyle 9&2&7\cr \frac{253}{134}&\frac{311}{67}&-\frac{503}{134}\cr -\frac{8080}{5253}&\frac{9494}{5253}&\frac{7676}{5253}} \right] \)，注意各向量彼此正交(垂直)。

各向量長度的平方為\( \matrix{\displaystyle \Vert\; b_1^*\Vert\;^2=9^2+2^2+7^2=134　　　　　　　　　　 \cr \Vert\; b_2^*\Vert\;^2=\left( \frac{253}{134} \right)^2+\left( \frac{311}{67} \right)^2+\left( -\frac{503}{134} \right)^2=\frac{5253}{134}　 \cr \Vert\; b_3^*\Vert\;^2=\left( -\frac{8080}{5253} \right)^2+\left( \frac{9494}{5253} \right)^2+\left( \frac{7676}{5253} \right)^2=\frac{40804}{5253}} \)
\( B_1=\left[ \matrix{9&2&7} \right] \)，\( b_1 \)向量長度就是\( \Vert\; b_1^* \Vert\; \)，平方為\( \Vert\; b_1^* \Vert\;^2=134 \)
\( B_2=\left[ \matrix{9&2&7 \cr 8&6&1} \right] \)，\( b_1,b_2 \)所張成的平行四邊形面積等同於\( b_1^*\)長、\(b_2^* \)寬所張成的長方形面積，平方為\( \displaystyle \Vert\; b_1^* \Vert\;^2 \cdot \Vert\; b_2^* \Vert\;^2=134 \cdot \frac{5253}{134}=5253 \)
\( B_3=\left[ \matrix{9&2&7 \cr 8&6&1 \cr 3&2&6} \right] \)，\( b_1,b_2,b_3 \)所張成的平行六面體體積等同於\( b_1^*\)長、\( b_2^*\)寬、\( b_3^* \)高所張成的長方體體積，平方為\( \displaystyle \Vert\; b_1^* \Vert\;^2 \cdot \Vert\; b_2^* \Vert\;^2 \cdot \Vert\; b_3^* \Vert\;^2=134 \cdot \frac{5253}{134} \cdot \frac{40804}{5253}=40804 \)
將全部sub lattice體積的平方全部乘起來\( 134 \cdot 5253 \cdot 40804=28722017208 \)就是Potential值


＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
以下程式計算Potential值

要先載入eigen.mac才能使用gramschmidt指令
(%i1)　load("eigen.mac");
(%o1)　C:/Program Files (x86)/Maxima-sbcl-5.37.1/share/maxima/5.37.1/share/matrix/eigen.mac

計算Potential值副程式
(%i2)
Pot(v):=block
(vstar:gramschmidt(v),
n:length(v[1]),
product((vstar[ i ].vstar[ i ])^(n-i+1),i,1,n)
)$

計算v的Potential值
(%i3)
v:matrix([9,2,7],
              [8,6,1],
              [3,2,6]);
Pot(v);
(%o3)　\( \left[ \matrix{9&2&7 \cr 8&6&1 \cr 3&2&6} \right] \)
(%o4)　28722017208



＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
以下程式為PotLLL執行結果

GramSchmit正交化
(%i1)
GramSchmit(v):=block
([m:length(v)],
vstar:copymatrix(v),
mu:zeromatrix(m,m),
for i:1 thru m do
  (for j:1 thru i do
     (mu[ i ][j]:v[ i ].vstar[j]/vstar[j].vstar[j],
      vstar[ i ]:v[ i ]-sum(mu[ i ][j]*vstar[j],j,1,i-1)
     )
  )
)$

檢查是否違反Size condition
(%i2)
SizeReduce(i,j):=block
([q:round(mu[ i ][j])],
if q#0 then
   (print("round(μ[",i,"][",j,"])=",round(mu[ i ][j]),"違反Size condition"),
    print("v=",v,"  第",i,"列-",q,"*第",j,"列  v=",rowop(v,i,j,q)),
    v:rowop(v,i,j,q)
   )
)$

第l列移動到第k列
(%i3)
MoveRow(v,k,l) := block
([n:length(v)],
head(v, n) :=if n > 0 then head(submatrix(length(v), v), n - 1) else v,
tail(v, n) :=if n > 0 then tail(submatrix(1, v), n - 1) else v,
v:addrow(head(v, n - k+1), v[l], tail(v, k-1)),
if k<l then v:submatrix(l+1,v) else v:submatrix(l,v)
)$

PotLLL主程式
(%i4)
PotLLL(v):=block(
[l:1,n:length(v),delta:3/4],
while l<=n do
  (if l>1 then
     (print("第",l,"列進行SizeReduce,μ=",mu)
     ),
   for i:1 thru l-1 do
     (SizeReduce(l,i)
     ),
   GramSchmit(v),
   P:1,
   Pmin:1,
   k:1,
   for j:l-1 thru 1 step -1 do
     (P: P*(vstar[l].vstar[l]+sum(mu[l,i]^2*vstar[ i ].vstar[ i ],i,j,l-1))/(vstar[j].vstar[j]),
      if P< Pmin then
        (k:j,
         Pmin: P
        )
     ),
   print("第",l,"列移到第",k,"列有最小的potential比值",Pmin),
   if delta> Pmin then
     (print("potential比值",Pmin,"比",delta,"小"),
      tempv:v,
      v:MoveRow(v,k,l),
      print(tempv,"第",l,"列移到第",k,"列",v),
      GramSchmit(v),
      print("l值",l,"=>",k),
      l:k
     )
   else
     (print("potential比值",Pmin,"比",delta,"大,l值",l,"=>",l+1),
      l:l+1
     )
  ),
return(v)
)$

(%i6)
v:matrix([9,2,7],
              [8,6,1],
              [3,2,6]);
PotLLL(v);
(%o5)　\( \left[ \matrix{9&2&7 \cr 8&6&1 \cr 3&2&6} \right] \)
第1列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值1=>2

第2列進行SizeReduce,μ=\( \left[ \matrix{\displaystyle 1&0&0 \cr \frac{91}{134}&1&0 \cr \frac{73}{134}&-\frac{1015}{5253}&1} \right] \)
round(μ[2][1])=1違反Size condition
\( v=\left[ \matrix{9&2&7 \cr 8&6&1 \cr 3&2&6} \right] \)  第2列\(-1*\)第1列  \( \left[ \matrix{9&2&7 \cr -1&4&-6 \cr 3&2&6} \right] \)

第2列移到第1列有最小的potential比值\( \displaystyle \frac{53}{134} \)
potential比值\( \displaystyle \frac{53}{134} \)比\( \displaystyle \frac{3}{4} \)小
\( \left[ \matrix{9&2&7\cr -1&4&-6\cr 3&2&6} \right] \)第2列移到第1列\( \left[ \matrix{-1&4&-6 \cr 9&2&7 \cr 3&2&6} \right] \)
l值2=>1

第1列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值1=>2

第2列進行SizeReduce,\( \mu=\left[ \matrix{\displaystyle 1&0&0\cr -\frac{43}{53}&1&0\cr -\frac{31}{53}&\frac{2536}{5253}&1} \right] \)
round(μ[2][1])=-1違反Size condition
\( v=\left[ \matrix{-1&4&-6\cr 9&2&7\cr 3&2&6} \right] \)  第2列\( - -1* \)第1列  \( v=\left[ \matrix{-1&4&-6\cr 8&6&1\cr 3&2&6} \right] \)

第2列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值2=>3

第3列進行SizeReduce,\( \mu=\left[ \matrix{\displaystyle 1&0&0\cr \frac{10}{53}&1&0\cr -\frac{31}{53}&\frac{2536}{5253}&1} \right] \)
round(μ[3][1])=-1違反Size condition
\( v=\left[ \matrix{-1&4&-6\cr 8&6&1\cr 3&2&6} \right] \)  第3列\( - -1* \)第1列  \( v=\left[ \matrix{-1&4&-6\cr 8&6&1\cr 2&6&0} \right] \)

第2列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值2=>3

第3列進行SizeReduce,\( \mu=\left[ \matrix{\displaystyle 1&0&0\cr -\frac{9}{20}&1&0\cr \frac{13}{10}&-\frac{186}{409}&1} \right] \)
round(μ[3][1])=1違反Size condition
\( v=\left[ \matrix{2&6&0\cr -3&-2&-6\cr 8&6&1} \right] \)  第3列\( -1* \)第1列  \( v=\left[ \matrix{2&6&0\cr -3&-2&-6\cr 6&0&1} \right] \)

第3列移到第1列有最小的potential比值\( \displaystyle \frac{65440}{278409} \)
potential 比值\( \displaystyle \frac{65440}{278409} \)比\( \displaystyle \frac{3}{4} \)小
\( \left[ \matrix{-1&4&-6 \cr 8&6&1 \cr 2&6&0} \right] \)第3列移到第1列\(  \left[ \matrix{2&6&0 \cr -1&4&-6 \cr 8&6&1} \right]\)
l值3=>1

第1列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值1=>2
第2列進行SizeReduce,\( \mu=\left[ \matrix{\displaystyle 1&0&0 \cr \frac{11}{20}&1&0 \cr \frac{13}{10}&-\frac{186}{409}&1} \right] \)
round(μ[2][1])=1違反Size condition
\( v=\left[ \matrix{2&6&0 \cr -1&4&-6 \cr 8&6&1} \right] \)  第2列\(-1*\)第1列  \( v=\left[ \matrix{2&6&0 \cr -3 & -2 &-6 \cr 8&6&1} \right] \)

第2列移到第1列有最小的potential比值1
potential比值1比\( \displaystyle \frac{3}{4} \)大,l值2=>3
第3列進行SizeReduce,\( \mu=\left[ \matrix{\displaystyle 1&0&0 \cr -\frac{9}{20}&1&0 \cr \frac{13}{10}&-\frac{186}{409}&1} \right] \)
round(μ[3][1])=1違反Size condition
\( v=\left[ \matrix{2&6&0 \cr -3&-2&-6 \cr 8&6&1} \right] \)  第3列\(-1*\)第1列  \( v=\left[ \matrix{2&6&0 \cr -3&-2&-6 \cr 6&0&1} \right] \)
第3列移到第1列有最小的potential比值\( \displaystyle \frac{6179}{8180} \)
potential比值\( \displaystyle \frac{6179}{8180} \)比\( \displaystyle \frac{3}{4} \)大,l值3=>4
(%o6)　\( \left[ \matrix{2&6&0 \cr -3&-2&-6 \cr 6&0&1} \right] \)


參考資料
PotLLL: A Polynomial Time Version of LLL With Deep Insertions
https://arxiv.org/abs/1307.7534
PDF檔https://arxiv.org/pdf/1307.7534v1.pdf

6-1.改進LLL方法－Jacobi方法

之前的Gaussian/Lagrange方法針對兩個向量進行化簡
https://math.pro/db/viewthread.php?tid=1876&page=1#pid10232
化簡後的向量符合(1)\(  \displaystyle \Bigg\vert\; \frac{v_1 \cdot v_2}{v_1 \cdot v_1} \Bigg\vert\;\le \frac{1}{2} \)　(2)\( \Vert\; v_1 \Vert\; \le \Vert\; v_2 \Vert\; \)兩個條件
擴展到更多向量就是將任兩個向量進行Gaussian/Lagrange化簡，任兩個向量滿足
(1)\(  \displaystyle \Bigg\vert\; \frac{v_i \cdot v_j}{v_i \cdot v_i} \Bigg\vert\;\le \frac{1}{2} \)　(2)\( \Vert\; v_i \Vert\; \le \Vert\; v_j \Vert\; \)，\( 1 \le i<j\le n \)兩個條件



(%i1)
Jacobi(v):=block
([i,j,n:length(v),t,End],
do
  (End:true,
   for i:1 thru n-1 do
     (for j:i+1 thru n do
        (if v[ i ].v[ i ]>v[j].v[j] then
           (print("第",i,"列長度",sqrt(v[ i ].v[ i ]),"比第",j,"列長度",sqrt(v[j].v[j]),"長"),
            print(v,"第",i,"列和第",j,"列交換",v:rowswap(v,i,j)),
            End:false
           ),
        t:round(v[ i ].v[j]/v[ i ].v[ i ]),
        if t#0 then
           (print("t=round(v",i,".v",j,"/v",i,".v",i,")=round(",v[ i ].v[j],"/",v[ i ].v[ i ],")=",t,"不等於0"),
            print(v,"第",j,"列-",t,"*第",i,"列",v:rowop(v,j,i,t)),
            End:false
           )
        )
     ),
   if End=true then
     (print("化簡後v=",v,
             ",",vi.vj/vi.vi,"=",genmatrix(lambda([i,j],if i<j then v[ i ].v[j]/v[ i ].v[ i ] else "*"),n,n),
             ",向量長度=",genmatrix(lambda([i,j],sqrt(v[ i ].v[ i ])),n,1)),
      print("所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束"),
      return(v)
     ),
   print("------------------------")
  )
)$

(%i2)
v:matrix([1,1,1],
             [-1,0,2],
             [3,5,6]);
(%o2)　\( \left[ \matrix{1&1&1\cr -1&0&2\cr 3&5&6} \right] \)

(%i3)　Jacobi(v);
\( \displaystyle t=round(v1.v3/v1.v1)=round(14/3)=5\)不等於0
\( \left[ \matrix{1&1&1\cr -1&0&2\cr 3&5&6} \right] \)第3列-5*第1列\( \left[ \matrix{1&1&1\cr -1&0&2\cr -2&0&1} \right] \)
\( \displaystyle t=round(v2.v3/v2.v2)=round(4/5)=1\)不等於0
\( \left[ \matrix{1&1&1\cr -1&0&2\cr -2&0&1} \right] \)第3列-1*第2列\( \left[ \matrix{1&1&1\cr -1&0&2\cr -1&0&-1} \right] \)
------------------------
第1列長度\( \sqrt{3}\)比第3列長度\( \sqrt{2}\)長
\( \left[ \matrix{1&1&1\cr -1&0&2\cr -1&0&-1} \right] \)第1列和第3列交換\( \left[ \matrix{-1&0&-1\cr -1&0&2\cr 1&1&1} \right] \)
\( \displaystyle t=round(v1.v3/v1.v1)=round(-2/2)=-1\)不等於0
\( \left[ \matrix{-1&0&-1\cr -1&0&2\cr 1&1&1} \right] \)第3列--1*第1列\( \left[ \matrix{-1&0&-1\cr -1&0&2\cr 0&1&0} \right] \)
第2列長度\( \sqrt{5}\)比第3列長度1長
\( \left[ \matrix{-1&0&-1\cr -1&0&2\cr 0&1&0} \right] \)第2列和第3列交換\( \left[ \matrix{-1&0&-1\cr 0&1&0\cr -1&0&2} \right] \)
------------------------
第1列長度\(\sqrt{2}\)比第2列長度1長
\( \left[ \matrix{-1&0&-1\cr 0&1&0\cr -1&0&2} \right] \)第1列和第2列交換\( \left[ \matrix{0&1&0\cr -1&0&-1\cr -1&0&2} \right] \)
------------------------
化簡後\( v=\left[ \matrix{0&1&0\cr -1&0&-1\cr -1&0&2} \right] \),\( \frac{vi.vj}{vi^2}=\left[ \matrix{*&0&0\cr *&*&-\frac{1}{2} \cr *&*&*} \right] \),向量長度\( =\left[ \matrix{1 \cr \sqrt{2} \cr \sqrt{5}} \right] \)
所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束

(%o3)　\( \left[ \matrix{0&1&0\cr -1&0&-1\cr -1&0&2} \right] \)

(%i4)
v:matrix([1,1,7,2],
         [9,8,4,6],
         [1,8,5,7],
         [2,3,1,1]);
(%o4)　\( \left[ \matrix{1&1&7&2\cr 9&8&4&6\cr 1&8&5&7\cr 2&3&1&1} \right] \)

(%i5)　Jacobi(v);
\( \displaystyle t=round(v1.v2/v1.v1)=round(57/55)=1\)不等於0
\( \left[ \matrix{1&1&7&2\cr 9&8&4&6\cr 1&8&5&7\cr 2&3&1&1}\right] \)第2列-1*第1列\( \left[ \matrix{1&1&7&2\cr 8&7&-3&4\cr 1&8&5&7\cr 2&3&1&1} \right] \)
\( \displaystyle t=round(v1.v3/v1.v1)=round(58/55)=1\)不等於0
\( \left[ \matrix{1&1&7&2\cr 8&7&-3&4\cr 1&8&5&7\cr 2&3&1&1}\right] \)第3列-1*第1列\( \left[ \matrix{1&1&7&2\cr 8&7&-3&4\cr 0&7&-2&5\cr 2&3&1&1} \right] \)
第1列長度\( \sqrt{55}\)比第4列長度\( \sqrt{15}\)長
\( \left[ \matrix{1&1&7&2\cr 8&7&-3&4\cr 0&7&-2&5\cr 2&3&1&1}\right] \)第1列和第4列交換\( \left[ \matrix{2&3&1&1\cr 8&7&-3&4\cr 0&7&-2&5\cr 1&1&7&2} \right] \)
\( \displaystyle t=round(v1.v4/v1.v1)=round(14/15)=1\)不等於0
\( \left[ \matrix{2&3&1&1\cr 8&7&-3&4\cr 0&7&-2&5\cr 1&1&7&2}\right] \)第4列-1*第1列\( \left[ \matrix{2&3&1&1\cr 8&7&-3&4\cr 0&7&-2&5\cr -1&-2&6&1} \right] \)
第2列長度\( \sqrt{138}\)比第3列長度\( \sqrt{78}\)長
\( \left[ \matrix{2&3&1&1\cr 8&7&-3&4\cr 0&7&-2&5\cr -1&-2&6&1}\right] \)第2列和第3列交換\( \left[ \matrix{2&3&1&1\cr 0&7&-2&5\cr 8&7&-3&4\cr -1&-2&6&1} \right] \)
\( \displaystyle t=round(v2.v3/v2.v2)=round(75/78)=1\)不等於0
\( \left[ \matrix{2&3&1&1\cr 0&7&-2&5\cr 8&7&-3&4\cr -1&-2&6&1}\right] \)第3列-1*第2列\( \left[ \matrix{2&3&1&1\cr 0&7&-2&5\cr 8&0&-1&-1\cr -1&-2&6&1} \right] \)
第2列長度\( \sqrt{78}\)比第4列長度\( \sqrt{42}\)長
\( \left[ \matrix{2&3&1&1\cr 0&7&-2&5\cr 8&0&-1&-1\cr -1&-2&6&1}\right] \)第2列和第4列交換\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 8&0&-1&-1\cr 0&7&-2&5} \right] \)
------------------------
\( \displaystyle t=round(v1.v3/v1.v1)=round(14/15)=1\)不等於0
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 8&0&-1&-1\cr 0&7&-2&5}\right] \)第3列-1*第1列\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 6&-3&-2&-2\cr 0&7&-2&5} \right] \)
\( \displaystyle t=round(v1.v4/v1.v1)=round(24/15)=2\)不等於0
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 6&-3&-2&-2\cr 0&7&-2&5}\right] \)第4列-2*第1列\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 6&-3&-2&-2\cr -4&1&-4&3} \right] \)
第3列長度\( \sqrt{53}\)比第4列長度\( \sqrt{42}\)長
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 6&-3&-2&-2\cr -4&1&-4&3}\right] \)第3列和第4列交換\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 6&-3&-2&-2} \right] \)
\( \displaystyle t=round(v3.v4/v3.v3)=round(-25/42)=-1\)不等於0
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 6&-3&-2&-2}\right] \)第4列--1*第3列\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 2&-2&-6&1} \right] \)
------------------------
\( \displaystyle t=round(v2.v4/v2.v2)=round(-33/42)=-1\)不等於0
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 2&-2&-6&1}\right] \)第4列--1*第2列\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 1&-4&0&2} \right] \)
第3列長度\( \sqrt{42}\)比第4列長度\( \sqrt{21}\)長
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr -4&1&-4&3\cr 1&-4&0&2}\right] \)第3列和第4列交換\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 1&-4&0&2\cr -4&1&-4&3} \right] \)
------------------------
\( \displaystyle t=round(v1.v3/v1.v1)=round(-8/15)=-1\)不等於0
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 1&-4&0&2\cr -4&1&-4&3}\right] \)第3列--1*第1列\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 3&-1&1&3\cr -4&1&-4&3} \right] \)
第2列長度\( \sqrt{42}\)比第3列長度\(2\sqrt{5}\)長
\( \left[ \matrix{2&3&1&1\cr -1&-2&6&1\cr 3&-1&1&3\cr -4&1&-4&3}\right] \)第2列和第3列交換\( \left[ \matrix{2&3&1&1\cr 3&-1&1&3\cr -1&-2&6&1\cr -4&1&-4&3} \right] \)
------------------------
化簡後\( v=\left[ \matrix{2&3&1&1\cr 3&-1&1&3\cr -1&-2&6&1\cr -4&1&-4&3} \right] \),\( \frac{vi.vj}{vi^2}=\left[ \matrix{ *&\frac{7}{15}&-\frac{1}{15}&-\frac{2}{5} \cr *&*&\frac{2}{5}&-\frac{2}{5} \cr *&*&*&-\frac{19}{42} \cr *&*&*&*} \right] \),向量長度\( =\left[ \matrix{\sqrt{15} \cr 2\sqrt{5} \cr \sqrt{42} \cr \sqrt{42}} \right] \)
所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束
(%o5)　\( \left[ \matrix{2&3&1&1\cr 3&-1&1&3\cr -1&-2&6&1\cr -4&1&-4&3} \right] \)

6-2.改進LLL方法－Jacobi方法

原本的Gauss/Lagrange方法

本論文的Gauss/Lagrange方法

1  \( do\{\; \) 2　\(if(v_1的長度>v_2的長度)\) 3　　\( \{\; v_1和v_2兩個向量交換 \}\; \) 4　\(t=\lfloor\; v_1.v_2/v_1.v_1 \rceil\; \) 5  \(if(t!=0)\) 6    \(\{\; v_2=v_2−tv_1; \}\;\) 7　 \(\}\;\) 8  \(while(t!=0);\)

1  \(G=V^TV\) 2  \(U=I_2\) 3  \( if \) \(G(1,1)<G(2,2)\) 4    \( swap \) \(G(:,1)\) \(and\) \(G(:,2)\) 5    \( swap \) \(G(1,: )\) \(and\) \(G(2,: )\) 6    \( swap \) \(Z(:,1)\) \(and\) \(G(:,2)\) 7  \(end\) 8 9  \(while\) \(G(1,1)>G(2,2)\) 10   \(q=\lfloor\; G(1,2)/G(2,2) \rceil\;\) 11   \(G(:,2)=G(:,2)-q \times G(:,1)\) 12   \(G(2,: )=G(2,: )-q \times G(1,: )\) 13   \(U(:,2)=U(:,2)-q \times U(:,1)\) 14 \(end\) 15 16 \(return\) \(U\)

註：說明論文所用的符號
1.\(G=V^TV\)代表gram矩陣，其中\(G(i,i)=v_i.v_i,G(i,j)=v_i.v_j\)，對角線元素為每個向量長度的平方，其他元素為兩向量的內積。
2.\(B\)代表lattice，為文章一致改用\(V\)表示lattice。
3.\(Z\)代表unimodular矩陣，但一般常用\(U\)表示unimodular矩陣，其中\(det(U)=+1或-1\)。

Jacobi方法以列計算

Jacobi方法以行計算

1  \(G=VV^T\) 2  \(U=I_n\) 3 4  \(while\) \(not\) \(all\) \(pairs\) \((v_i,v_j)\) \(satisfy\) \( \Vert\; v_i \Vert\; \le \Vert\; v_j \Vert\; \) and \( \displaystyle |\; v_i v_j |\; \le \frac{\Vert\; v_i \Vert\;^2}{2} \) 5    \(for\) \(i=1\) \(to\) \(n-1\) 6      \(for\) \(j=i+1\) \(to\) \(n\) 7        \(q=G(i,j)/G(i,i)\) 8        \(if\) \(|\; q |\;>1/2\) 9          \(G(:,j)=G(:,j)-\lfloor\; q \rceil\; \times G(:,i)\) 10         \(G(j,: )=G(j,: )-\lfloor\; q \rceil\; \times G(i,: )\) 11         \(U(j,: )=U(j,: )-\lfloor\; q \rceil\; \times U(i,: )\) 12       \(end\) 13       \(if\) \(G(i,i)>G(j,j)\) 14         \(swap\) \(G(:,i)\) \(and\) \(G(:,j)\) 15         \(swap\) \(G(i,: )\) \(and\) \(G(j,: )\) 16         \(swap\) \(U(i,: )\) \(and\) \(U(j,: )\) 17       \(end\) 18     \(end\) 19   \(end\) 20 \(end\) 21 22 \(return\) \(UV\)

1  \(G=V^TV\) 2  \(U=I_n\) 3 4  \(while\) \(not\) \(all\) \(pairs\) \((v_i,v_j)\) \(satisfy\) \( \Vert\; v_i \Vert\; \le \Vert\; v_j \Vert\; \) and \( \displaystyle |\; v_i v_j |\; \le \frac{\Vert\; v_i \Vert\;^2}{2} \) 5    \(for\) \(i=1\) \(to\) \(n-1\) 6      \(for\) \(j=i+1\) \(to\) \(n\) 7        \(q=G(i,j)/G(i,i)\) 8        \(if\) \(|\; q |\;>1/2\) 9          \(G(:,j)=G(:,j)-\lfloor\; q \rceil\; \times G(:,i)\) 10         \(G(j,: )=G(j,: )-\lfloor\; q \rceil\; \times G(i,: )\) 11         \(U(:,j)=U(:,j)-\lfloor\; q \rceil\; \times U(:,i)\) 12       \(end\) 13       \(if\) \(G(i,i)>G(j,j)\) 14         \(swap\) \(G(:,i)\) \(and\) \(G(:,j)\) 15         \(swap\) \(G(i,: )\) \(and\) \(G(j,: )\) 16         \(swap\) \(U(:,i)\) \(and\) \(U(:,j)\) 17       \(end\) 18     \(end\) 19   \(end\) 20 \(end\) 21 22 \(return\) \(VU\)

註：
論文演算法以行向量運算，本文章卻以列向量運算，所以將兩個演算法並列出來，其中第1,11,16,22行有差異。
論文演算法第7行\(q=G(i,j)/G(j,j)\)會導致lattice尚未化簡完程式就提早結束，更正為\(q=G(i,j)/G(i,i)\)。

參考資料：
http://www.cas.mcmaster.ca/~qiao/
Filip Jeremic and Sanzheng Qiao. A Parallel Jacobi-type Lattice Basis Reduction Algorithm. International Journal of Numerical Analysis and Modeling, Series B. Vol. 5, No. 1-2, 2014, 1-12.
http://www.cas.mcmaster.ca/~qiao/publications/JQ14.pdf
Filip Jeremic and Sanzheng Qiao. A GPU Implementation of a Jacobi Method for Lattice Basis Reduction. Proceedings of International Workshop on Data-Intensive Scientific Discovery (DISD), 2013, Shanghai University, Shanghai, China, August 1-4, 2013.
http://www.cas.mcmaster.ca/~qiao/publications/JQ13.pdf



以列向量運算的Jacobi方法
(%i1)
JacobibyRow(V):=block
([End,G,n:length(V),U,i,j],
G:V.transpose(V),
U:ident(n),
do
(End:true,
for i:1 thru n-1 do
  (for j:i+1 thru n do
     (q:G[i,j]/G[i,i],
      if abs(q)>1/2 then
        (print("G(",i,",",j,")/G(",i,",",i,")=",q,">",1/2,"  ,  round(",q,")=",round(q)),
         print("G=",G,"  第",j,"行-",round(q),"*第",i,"行  G=",G:columnop(G,j,i,round(q))),
         print("G=",G,"  第",j,"列-",round(q),"*第",i,"列  G=",G:rowop(G,j,i,round(q))),
         print("U=",U,"  第",j,"列-",round(q),"*第",i,"列  U=",U:rowop(U,j,i,round(q))),
         End:false
        ),
      if G[i,i]>G[j,j] then
        (print("G(",i,",",i,")=",G[i,i],">G(",j,",",j,")=",G[j,j]),
         print("G=",G,"  第",i,"行和第",j,"行交換  G=",G:columnswap(G,i,j)),
         print("G=",G,"  第",i,"列和第",j,"列交換  G=",G:rowswap(G,i,j)),
         print("U=",U,"  第",i,"列和第",j,"列交換  U=",U:rowswap(U,i,j)),
         End:false
        )
     )
    ),
   if End=true then
     (print("G=",G,",",vi.vj/vi.vi,"=",genmatrix(lambda([i,j],if i<j then G[i,j]/G[i,i] else "*"),n,n),
             ",向量長度=",genmatrix(lambda([i,j],sqrt(G[i,i])),n,1)),
      print("所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束"),
      print("化簡後向量V=U.V=",U,V,"=",U.V),
      return(U.V)
     ),
    print("------------------------")
  )
)$

列向量的lattice
(%i2)
V:matrix([1,1,1],
             [-1,0,2],
             [3,5,6]);

(%o2)　\( \left[ \matrix{1&1&1\cr -1&0&2\cr 3&5&6} \right] \)

(%i3)　JacobibyRow(V);
\( \displaystyle G( 1,3)/G(1,1)=\frac{14}{3}>\frac{1}{2}\)  ,  \( \displaystyle round(\frac{14}{3})=5\)
\( G=\left[ \matrix{3&1&14\cr 1&5&9\cr 14&9&70} \right] \)第3行-5*第1行  \( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr 14&9&0} \right] \)
\( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr 14&9&0} \right] \)第3列-5*第1列  \( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr -1&4&5} \right] \)
\( U=\left[ \matrix{1&0&0\cr 0&1&0\cr 0&0&1} \right] \)第3列-5*第1列  \( U=\left[ \matrix{1&0&0\cr 0&1&0\cr -5&0&1} \right] \)
\( \displaystyle G( 2,3)/G(2,2)=\frac{4}{5}>\frac{1}{2}\)  ,  \( \displaystyle round(\frac{4}{5})=1\)
\( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr -1&4&5} \right] \)第3行-1*第2行  \( G=\left[ \matrix{3&1&-2\cr
1&5&-1\cr -1&4&1} \right] \)
\( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -1&4&1} \right] \)第3列-1*第2列  \( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -2&-1&2} \right] \)
\( U=\left[ \matrix{1&0&0\cr 0&1&0\cr -5&0&1} \right] \)第3列-1*第2列  \( U=\left[ \matrix{1&0&0\cr
0&1&0\cr -5&-1&1} \right] \)
\( G( 2,2)=5>G(3,3)=2 \)
\( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -2&-1&2} \right] \)第2行和第3行交換  \( G=\left[ \matrix{3&-2&1\cr 1&-1&5\cr -2&2&-1} \right] \)
\( G=\left[ \matrix{3&-2&1\cr 1&-1&5\cr -2&2&-1} \right] \)第2列和第3列交換  \( G=\left[ \matrix{3&-2&1\cr -2&2&-1\cr 1&-1&5} \right] \)
\( U=\left[ \matrix{1&0&0\cr 0&1&0\cr -5&-1&1} \right] \)第2列和第3列交換  \( U=\left[ \matrix{1&0&0\cr -5&-1&1\cr 0&1&0} \right] \)
------------------------
\( \displaystyle G( 1,2)/G(1,1)=-\frac{2}{3}>\frac{1}{2}\)  ,  \( \displaystyle round(\frac{-2}{3})=-1\)
\( G=\left[ \matrix{3&-2&1\cr -2&2&-1\cr 1&-1&5} \right] \)第2行--1*第1行  \( G=\left[ \matrix{3&1&1\cr -2&0&-1\cr 1&0&5} \right] \)
\( G=\left[ \matrix{3&1&1\cr -2&0&-1\cr 1&0&5} \right] \)第2列--1*第1列  \( G=\left[ \matrix{3&1&1\cr
1&1&0\cr 1&0&5} \right] \)
\( U=\left[ \matrix{1&0&0\cr -5&-1&1\cr 0&1&0} \right] \)第2列--1*第1列  \( U=\left[ \matrix{1&0&0\cr -4&-1&1\cr 0&1&0} \right] \)
\( G( 1,1)=3>G(2,2)=1\)
\( G=\left[ \matrix{3&1&1\cr 1&1&0\cr 1&0&5} \right] \)第1行和第2行交換  \( G=\left[ \matrix{1&3&1\cr 1&1&0\cr 0&1&5} \right] \)
\( G=\left[ \matrix{1&3&1\cr 1&1&0\cr 0&1&5} \right] \)第1列和第2列交換  \( G=\left[ \matrix{1&1&0\cr 1&3&1\cr 0&1&5} \right] \)
\( U=\left[ \matrix{1&0&0\cr -4&-1&1\cr 0&1&0} \right] \)第1列和第2列交換  \( U=\left[ \matrix{
-4&-1&1\cr 1&0&0\cr 0&1&0} \right] \)
------------------------
\( \displaystyle G( 1,2)/G(1,1)=1>\frac{1}{2}\)  ,  \( \displaystyle round(1)=1\)
\( G=\left[ \matrix{1&1&0\cr 1&3&1\cr 0&1&5} \right] \)第2行-1*第1行  \( G=\left[ \matrix{1&0&0\cr
1&2&1\cr 0&1&5} \right] \)
\( G=\left[ \matrix{1&0&0\cr 1&2&1\cr 0&1&5} \right] \)第2列-1*第1列  \( G=\left[ \matrix{1&0&0\cr
0&2&1\cr 0&1&5} \right] \)
\( U=\left[ \matrix{-4&-1&1\cr 1&0&0\cr 0&1&0} \right] \)第2列-1*第1列  \( U=\left[ \matrix{-4&-1&1\cr 5&1&-1\cr 0&1&0} \right] \)
------------------------
\( G=\left[ \matrix{1&0&0\cr 0&2&1\cr 0&1&5} \right]\),\( \displaystyle \frac{vi.vj}{vi^2}=\left[\matrix{\displaystyle *&0&0\cr *&*&\frac{1}{2}\cr *&*&*}\right]\),向量長度\(=\left[ \matrix{1\cr \sqrt{2}\cr \sqrt{5}} \right] \)
所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束
化簡後向量\( V=U.V=\left[ \matrix{-4&-1&1\cr 5&1&-1\cr 0&1&0} \right] \left[\matrix{1&1&1\cr -1&0&2\cr 3&5&6}\right]=\left[ \matrix{0&1&0\cr 1&0&1\cr -1&0&2} \right]\)
(%o3)　\( \left[ \matrix{0&1&0\cr 1&0&1\cr -1&0&2} \right] \)



－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－



以行向量運算的Jacobi方法
(%i1)
JacobibyColumn(V):=block
([End,G,n:length(V),U,i,j],
G:transpose(V).V,
U:ident(n),
do
(End:true,
for i:1 thru n-1 do
  (for j:i+1 thru n do
     (q:G[i,j]/G[i,i],
      if abs(q)>1/2 then
        (print("G(",i,",",j,")/G(",i,",",i,")=",q,">",1/2,"  ,  round(",q,")=",round(q)),
         print("G=",G,"  第",j,"行-",round(q),"*第",i,"行  G=",G:columnop(G,j,i,round(q))),
         print("G=",G,"  第",j,"列-",round(q),"*第",i,"列  G=",G:rowop(G,j,i,round(q))),
         print("U=",U,"  第",j,"行-",round(q),"*第",i,"行  U=",U:columnop(U,j,i,round(q))),
         End:false
        ),
      if G[i,i]>G[j,j] then
        (print("G(",i,",",i,")=",G[i,i],">G(",j,",",j,")=",G[j,j]),
         print("G=",G,"  第",i,"行和第",j,"行交換  G=",G:columnswap(G,i,j)),
         print("G=",G,"  第",i,"列和第",j,"列交換  G=",G:rowswap(G,i,j)),
         print("U=",U,"  第",i,"行和第",j,"行交換  U=",U:columnswap(U,i,j)),
         End:false
        )
     )
    ),
   if End=true then
     (print("G=",G,",",vi.vj/vi.vi,"=",genmatrix(lambda([i,j],if i<j then G[i,j]/G[i,i] else "*"),n,n),
             ",向量長度=",genmatrix(lambda([i,j],sqrt(G[i,i])),n,1)),
      print("所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束"),
      print("化簡後向量V=V.U=",V,U,"=",V.U),
      return(V.U)
     ),
    print("------------------------")
  )
)$

行向量的lattice
(%i2)
V:matrix([1,-1,3],
             [1,0,5],
             [1,2,6]);
(%o2)　\( \left[ \matrix{1&-1&3\cr 1&0&5\cr 1&2&6} \right] \)

(%i3)　JacobibyColumn(V);
\( \displaystyle G(1,3)/G(1,1)=\frac{14}{3}>\frac{1}{2} \)  ,  \( \displaystyle round(14/3)=5 \)
\( G=\left[ \matrix{3&1&14\cr 1&5&9\cr 14&9&70} \right] \)  第3行-5*第1行  \( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr 14&9&0} \right] \)
\( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr 14&9&0} \right] \)  第3列-5*第1列  \( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr -1&4&5} \right] \)
\( U=\left[ \matrix{1&0&0\cr 0&1&0\cr 0&0&1} \right] \)  第3行-5*第1行  \( U=\left[ \matrix{1&0&-5\cr 0&1&0\cr 0&0&1} \right] \)
\( \displaystyle G(2,3)/G(2,2)=\frac{4}{5}>\frac{1}{2} \)  ,  \( \displaystyle round(4/5)=1\)
\( G=\left[ \matrix{3&1&-1\cr 1&5&4\cr -1&4&5} \right] \)  第3行-1*第2行  \( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -1&4&1} \right] \)
\( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -1&4&1} \right] \)  第3列-1*第2列  \( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -2&-1&2} \right] \)
\( U=\left[ \matrix{1&0&-5\cr 0&1&0\cr 0&0&1} \right] \)  第3行-1*第2行  \( U=\left[ \matrix{1&0&-5\cr 0&1&-1\cr 0&0&1} \right] \)
\( G(2,2)=5>G(3,3)=2\)
\( G=\left[ \matrix{3&1&-2\cr 1&5&-1\cr -2&-1&2} \right] \)  第2行和第3行交換  \( G=\left[ \matrix{3&-2&1\cr 1&-1&5\cr -2&2&-1} \right] \)
\( G=\left[ \matrix{3&-2&1\cr 1&-1&5\cr -2&2&-1} \right] \)  第2列和第3列交換  \( G=\left[ \matrix{3&-2&1\cr -2&2&-1\cr 1&-1&5} \right] \)
\( U=\left[ \matrix{1&0&-5\cr 0&1&-1\cr 0&0&1} \right] \)  第2行和第3行交換  \( U=\left[ \matrix{
1&-5&0\cr 0&-1&1\cr 0&1&0} \right] \)
------------------------
\( \displaystyle G(1,2)/G(1,1)=-\frac{2}{3}>\frac{1}{2} \)  ,  \( \displaystyle round(-\frac{2}{3})=-1\)
\( G=\left[ \matrix{3&-2&1\cr -2&2&-1\cr 1&-1&5} \right] \)  第2行--1*第1行  \( G=\left[ \matrix{3&1&1\cr -2&0&-1\cr 1&0&5} \right] \)
\( G=\left[ \matrix{3&1&1\cr -2&0&-1\cr 1&0&5} \right] \)  第2列--1*第1列  \( G=\left[ \matrix{3&1&1\cr 1&1&0\cr 1&0&5} \right] \)
\( U=\left[ \matrix{1&-5&0\cr 0&-1&1\cr 0&1&0} \right] \)  第2行--1*第1行  \( U=\left[ \matrix{1&-4&0\cr 0&-1&1\cr 0&1&0} \right] \)
\( G(1,1)=3>G(2,2)=1\)
\( G=\left[ \matrix{3&1&1\cr 1&1&0\cr 1&0&5} \right] \)  第1行和第2行交換  \( G=\left[ \matrix{1&3&1\cr 1&1&0\cr 0&1&5} \right] \)
\( G=\left[ \matrix{1&3&1\cr 1&1&0\cr 0&1&5} \right] \)  第1列和第2列交換  \( G=\left[ \matrix{
1&1&0\cr 1&3&1\cr 0&1&5} \right] \)
\( U=\left[ \matrix{1&-4&0\cr 0&-1&1\cr 0&1&0} \right] \)  第1行和第2行交換  \( U=\left[ \matrix{
-4&1&0\cr -1&0&1\cr 1&0&0} \right] \)
------------------------
\( \displaystyle G(1,2)/G(1,1)=1>\frac{1}{2} \)  ,  \( round(1)=1\)
\( G=\left[ \matrix{1&1&0\cr 1&3&1\cr 0&1&5} \right] \)  第2行-1*第1行  \( G=\left[ \matrix{1&0&0\cr 1&2&1\cr 0&1&5} \right] \)
\( G=\left[ \matrix{1&0&0\cr 1&2&1\cr 0&1&5} \right] \)  第2列-1*第1列  \( G=\left[ \matrix{1&0&0\cr 0&2&1\cr 0&1&5} \right] \)
\( U=\left[ \matrix{-4&1&0\cr -1&0&1\cr 1&0&0} \right] \)  第2行-1*第1行  \( U=\left[ \matrix{-4&5&0\cr -1&1&1\cr 1&-1&0} \right] \)
------------------------
\( G=\left[ \matrix{1&0&0\cr 0&2&1\cr 0&1&5}\right],\frac{vi.vj}{vi^2}=\left[ \matrix{*&0&0\cr *&*&\frac{1}{2}\cr *&*&*}\right] \),向量長度\( =\left[ \matrix{1\cr \sqrt{2}\cr \sqrt{5}} \right] \)
所有向量符合(1)|vi.vj/vi.vi|≦1/2，(2)∥vi∥≦∥vj∥，1≦i<j≦n兩個條件，程式結束
化簡後向量\( V=V.U=\left[ \matrix{1&-1&3\cr 1&0&5\cr 1&2&6} \right] \left[ \matrix{-4&5&0\cr -1&1&1\cr 1&-1&0} \right]=\left[ \matrix{0&1&-1\cr 1&0&0\cr 0&1&2} \right]\)
(%o3)　\( \left[ \matrix{0&1&-1\cr 1&0&0\cr 0&1&2} \right] \)

[ 本帖最後由 bugmens 於 2017-12-31 13:54 編輯 ]