色々と修正

Author: satoshin-des

PLM.py

@@ -10,7 +10,6 @@
 """
 
 import numpy as np
-import sympy
 import sys
 
 class lattice():
@@ -56,6 +55,19 @@ def volume(self) -> float:
         return np.sqrt(np.linalg.det(self.basis * self.basis.T))
 
 
+    def dual(self):
+        """Computes dual-lattice.
+
+        Returns:
+            lattice: Dual-lattice.
+        """
+        c = lattice(self.basis)
+        self.basis_star, self.mu = self.GSO()
+        c.basis = np.matmul(np.matmul(np.linalg.inv(self.mu.T), np.linalg.inv(np.matmul(self.basis_star, self.basis_star.T))), self.basis_star)
+        c.basis_star, c.mu = c.GSO()
+        return c
+
+
     def GSO(self, mode: str = "normal") -> np.ndarray:
         """Gram-Schmidt's method.
 
@@ -82,7 +94,6 @@ def GSO(self, mode: str = "normal") -> np.ndarray:
                 self.basis_star[i] -= self.mu[i, j] * np.copy(self.basis_star[j])
             if mode == "square" or mode == "both":
                 self.B[i] = np.dot(self.basis_star[i], self.basis_star[i])
-       
         if mode == "normal":
             return self.basis_star, self.mu
         elif mode == "square":
@@ -171,7 +182,7 @@ def Gauss(self) -> np.ndarray:
 
 
     def Lagrange(self) -> np.ndarray:
-        """Lagrange reduction
+        """Lagrange reduction. This is other name of Gauss().
 
         Returns:
             np.ndarray: Lagrange reduced lattice basis matrix.
@@ -199,15 +210,13 @@ def LLL(self, delta: float = 0.99) -> np.ndarray:
             np.ndarray: LLL-reduced basis
         """
         self.B, self.mu = self.GSO(mode = "square")
-
         k = 1
         while k < self.nrows:
             for j in range(k)[::-1]:
                 if abs(self.mu[k, j]) > 0.5:
                     q = round(self.mu[k, j])
                     self.basis[k] -= q * np.copy(self.basis[j])
                     self.mu[k, : j + 1] -= q * np.copy(self.mu[j, : j + 1])
-
             if self.B[k] >= (delta - self.mu[k, k - 1] * self.mu[k, k - 1]) * self.B[k - 1]:
                 k += 1
             else:
@@ -257,8 +266,19 @@ def DeepLLL(self, delta : float = 0.99) -> np.ndarray:
                     self.B, self.mu = self.GSO(mode = "square")
                     k = max(i - 1, 0)
             k += 1
-        
         return self.basis
+    
+
+    def DEEP(self, delta: float = 0.99) -> np.ndarray:
+        """Other name of DeepLLL().
+
+        Args:
+            delta (float, optional): Reduction parameter. Defaults to 0.99.
+
+        Returns:
+            np.ndarray: Deep-LLL-reduced lattice basis matrix.
+        """
+        return self.DeepLLL(delta = delta)
 
 
     def PotLLL(self, delta: float = 0.99) -> np.ndarray:
@@ -281,8 +301,7 @@ def PotLLL(self, delta: float = 0.99) -> np.ndarray:
                     self.mu[l, : j + 1] -= q * np.copy(self.mu[j, : j + 1])
             P = P_min = 1.; k = 0
             for j in range(l)[::-1]:
-                S = np.sum(self.mu[l, j: l] * self.mu[l, j: l] * self.B[j: l])
-                P *= (self.B[l] + S) / self.B[j]
+                P *= (self.B[l] + np.sum(self.mu[l, j: l] * self.mu[l, j: l] * self.B[j: l])) / self.B[j]
                 if P < P_min: k = j; P_min = P
             if delta > P_min:
                 t = np.copy(self.basis[l])
@@ -370,7 +389,7 @@ def _ENUM_(mu: np.ndarray, B: np.ndarray, n: int, delta: float) -> np.ndarray:
                     w[k] = 1
                 else:
                     k += 1
-                    if k == n: return np.zeros(n, int)
+                    if k == n: return None
                     r[k - 1] = k
                     if k >= last_nonzero:
                         last_nonzero = k
@@ -386,7 +405,7 @@ def _ENUM_(mu: np.ndarray, B: np.ndarray, n: int, delta: float) -> np.ndarray:
         while True:
             pre_ENUM_v = np.copy(ENUM_v)
             ENUM_v = _ENUM_(self.mu, self.B, self.nrows, delta)
-            if np.all(ENUM_v == 0): return pre_ENUM_v
+            if ENUM_v is None: return np.matmul(pre_ENUM_v, self.basis)
             delta *= 0.99
 
 
@@ -417,45 +436,6 @@ def project_basis(self, k: int, l: int) -> np.ndarray:
         return pi_b
 
 
-    def BKZ(self, beta: int, delta: float = 0.99) -> np.ndarray:
-        """BKZ-reduces a lattice basis matrix(algorithm is from C. P. Schnorr and M. Euchner(1994)).
-
-        Now this function is very unstable.
-
-        Args:
-            beta (int): Block size.
-            delta (float, optional): Reduction parameter. Defaults to 0.99.
-
-        Returns:
-            np.ndarray: BKZ-reduced basis matirx.
-        """
-        self.basis = self.LLL(delta = delta)
-        z = k = 0
-        while z < self.nrows - 2:
-            print(z)
-            if k == self.nrows - 1: k = 0
-            k1 = k; k += 1
-            l = min(k1 + beta, self.nrows); h = min(l + 1, self.nrows)
-            self.B, self.mu = self.GSO(mode = "square")
-            p = lattice(self.basis[k1: l, k1: l])
-            w = p.ENUM_SVP(); s = w @ self.project_basis(k1, l - 1)
-            if (not np.all(s == 0)) and self.B[k1] > np.dot(s, s):
-                z = 0
-                c = lattice(np.zeros((h + 1, self.ncols)))
-                c.basis[: k1] = np.copy(self.basis[: k1])
-                c.basis[k1] = w @ self.basis[k1: l]
-                c.basis[k: h + 1] = np.copy(self.basis[k1: h])
-                _, inds = sympy.Matrix(c.basis).T.rref()
-                c.basis = np.copy(c.basis[np.array(inds)]); c.nrows = h
-                c.basis = c.LLL(delta = delta)
-                self.basis[: h] = np.copy(c.basis)
-            else:
-                z += 1
-                c = lattice(self.basis[: h])
-                self.basis[: h] = c.LLL(delta = delta)
-        return self.basis
-
-
     def Babai(self, w: np.ndarray):
         """Computes an approximate solution of CVP for target w using Babai's nearest plane algorithm(algorithm is from L. Babai(1986)).
 
@@ -519,7 +499,7 @@ def _ENUM_(mu: np.ndarray, B: np.ndarray, n: int, a: np.ndarray, R: float) -> np
                     w[k] = 1
                 else:
                     k += 1
-                    if k == n: return np.zeros(n, int)
+                    if k == n: return None
                     r[k - 1] = k
                     if v[k] > c[k]:	v[k] -= w[k]
                     else: v[k] += w[k]
@@ -528,11 +508,11 @@ def _ENUM_(mu: np.ndarray, B: np.ndarray, n: int, a: np.ndarray, R: float) -> np
         self.B, self.mu = self.GSO(mode = "square")
         ENUM_v = np.zeros(self.nrows, int)
         babaivec = self.Babai(t); R = np.dot(babaivec - t, babaivec - t)
-        a = t @ np.linalg.pinv(self.basis)
+        a = np.matmul(t, np.linalg.pinv(self.basis))
         while True:
             pre_ENUM_v = np.copy(ENUM_v)
             ENUM_v = _ENUM_(self.mu, self.B, self.nrows, a, R)
-            if np.all(ENUM_v == 0): return pre_ENUM_v @ self.basis
+            if ENUM_v is None: return np.matmul(pre_ENUM_v, self.basis)
             R *= 0.99
 
 

main.py

@@ -3,14 +3,14 @@
 import sympy
 
 def main():
-    b = PLM.random_lattice(40)
+    b = PLM.random_lattice(70)
     print(b.basis)
-    a = b.PotLLL()[0]; print(a, np.linalg.norm(a))
+    #a = b.PotLLL()[0]; print(a, np.linalg.norm(a))
     #a = PLM.lattice.ENUM(b); print(a @ b.basis, np.linalg.norm(a @ b.basis))
-    #a = b.BKZ(15); print(a, np.linalg.norm(a[0]))
-    t = np.random.randint(10, 20, size=40)
+    t = np.random.randint(10, 20, size=70)
     print(t)
     print(b.ENUM_CVP(t))
+    print(b.ENUM_SVP())
 
 
 if __name__ == '__main__':