diff --git a/python/examples/repetition_code.py b/python/examples/repetition_code.py index 812aa87..5c64be1 100644 --- a/python/examples/repetition_code.py +++ b/python/examples/repetition_code.py @@ -69,51 +69,59 @@ from hccd import path_tracker, homotopy_generator def track_path(args): - H = np.array([[1, 1, 1]]) + H = np.array([[1, 1, 0, 0], + [0, 1, 1, 0], + [0, 0, 1, 1]]) - # TODO: Make this work homotopy = homotopy_generator.HomotopyGenerator(H) + print(f"G: {homotopy.G}") + print(f"F: {homotopy.F}") + print(f"H: {homotopy.H}") + print(f"DH: {homotopy.DH}") + tracker = path_tracker.PathTracker(homotopy, args.euler_step_size, args.euler_max_tries, args.newton_max_iter, args.newton_convergence_threshold, args.sigma) ys_start, ys_prime, ys_hat_e, ys = [], [], [], [] - try: - y = np.zeros(3) - for i in range(args.num_iterations): - y_start, y_prime, y_hat_e, y = tracker.transparent_step(y) - ys_start.append(y_start) - ys_prime.append(y_prime) - ys_hat_e.append(y_hat_e) - ys.append(y) - print(f"Iteration {i}: {y}") - except Exception as e: - print(f"Error: {e}") + y = np.zeros(5) + np.array([0.5, 0.48, -0.1, 0.2, 0 ]) + for i in range(args.num_iterations): + y_start, y_prime, y_hat_e, y = tracker.transparent_step(y) + ys_start.append(y_start) + ys_prime.append(y_prime) + ys_hat_e.append(y_hat_e) + ys.append(y) + print(f"Iteration {i}: {y}") ys_start = np.array(ys_start) ys_prime = np.array(ys_prime) ys_hat_e = np.array(ys_hat_e) ys = np.array(ys) - df = pd.DataFrame({"x1b": ys_start[:, 0], - "x2b": ys_start[:, 1], - "tb": ys_start[:, 2], - "x1p": ys_prime[:, 0], - "x2p": ys_prime[:, 1], - "tp": ys_prime[:, 2], - "x1e": ys_hat_e[:, 0], - "x2e": ys_hat_e[:, 1], - "te": ys_hat_e[:, 2], - "x1n": ys[:, 0], - "x2n": ys[:, 1], - "tn": ys[:, 2] - }) + for y in np.transpose(ys): + plt.plot(y) - if args.output: - df.to_csv(args.output, index=False) - else: - print(df) + plt.show() + + # df = pd.DataFrame({"x1b": ys_start[:, 0], + # "x2b": ys_start[:, 1], + # "tb": ys_start[:, 2], + # "x1p": ys_prime[:, 0], + # "x2p": ys_prime[:, 1], + # "tp": ys_prime[:, 2], + # "x1e": ys_hat_e[:, 0], + # "x2e": ys_hat_e[:, 1], + # "te": ys_hat_e[:, 2], + # "x1n": ys[:, 0], + # "x2n": ys[:, 1], + # "tn": ys[:, 2] + # }) + # + # if args.output: + # df.to_csv(args.output, index=False) + # else: + # print(df) def main(): diff --git a/python/hccd/homotopy_generator.py b/python/hccd/homotopy_generator.py index 79f6264..bd6f423 100644 --- a/python/hccd/homotopy_generator.py +++ b/python/hccd/homotopy_generator.py @@ -17,78 +17,63 @@ class HomotopyGenerator: self.H_matrix = parity_check_matrix self.num_checks, self.num_vars = parity_check_matrix.shape - # Create symbolic variables self.x_vars = [sp.symbols(f'x{i+1}') for i in range(self.num_vars)] self.t = sp.symbols('t') - # Generate G, F, and H self.G = self._create_G() self.F = self._create_F() self.H = self._create_H() + self.DH = self._create_DH(self.H) - # Convert to callable functions self._H_lambda = self._create_H_lambda() self._DH_lambda = self._create_DH_lambda() def _create_G(self) -> List[sp.Expr]: - """Create G polynomial system (the starting system).""" - # For each variable xi, add the polynomial [xi] G = [] for var in self.x_vars: G.append(var) return G - def _create_F(self) -> List[sp.Expr]: - """Create F polynomial system (the target system).""" + def _create_F(self) -> sp.MutableMatrix: F = [] - # Add 1 - xi^2 for each variable for var in self.x_vars: F.append(1 - var**2) - # Add parity check polynomials: 1 - x1*x2*...*xk for each parity check for row in self.H_matrix: - # Create product of variables that participate in this check term = 1 for i, bit in enumerate(row): if bit == 1: term *= self.x_vars[i] - if term != 1: # Only add if there are variables in this check - F.append(1 - term) + F.append(1 - term) - return F + groebner_basis = sp.groebner(F) + return sp.MutableMatrix(groebner_basis) def _create_H(self) -> List[sp.Expr]: - """Create the homotopy H = (1-t)*G + t*F.""" H = [] - # Make sure G and F have the same length - # Repeat variables from G if needed to match F's length - G_extended = self.G.copy() - while len(G_extended) < len(self.F): - # Cycle through variables to repeat - for i in range(min(self.num_vars, len(self.F) - len(G_extended))): - G_extended.append(self.x_vars[i % self.num_vars]) - - # Create the homotopy - for g, f in zip(G_extended, self.F): + for g, f in zip(self.G, self.F): H.append((1 - self.t) * g + self.t * f) return H + def _create_DH(self, H: List[sp.Expr]) -> sp.MutableMatrix: + all_vars = self.x_vars + [self.t] + DH = sp.Matrix([[sp.diff(expr, var) + for var in all_vars] for expr in self.H]) + + return DH + def _create_H_lambda(self) -> Callable: - """Create a lambda function to evaluate H.""" all_vars = self.x_vars + [self.t] return sp.lambdify(all_vars, self.H, 'numpy') def _create_DH_lambda(self) -> Callable: - """Create a lambda function to evaluate the Jacobian of H.""" all_vars = self.x_vars + [self.t] - jacobian = sp.Matrix([[sp.diff(expr, var) - for var in all_vars] for expr in self.H]) - return sp.lambdify(all_vars, jacobian, 'numpy') + return sp.lambdify(all_vars, self.DH, 'numpy') def evaluate_H(self, y: np.ndarray) -> np.ndarray: """