From ec598166e5bfcb15cb124bf3e62904f42596f5c6 Mon Sep 17 00:00:00 2001 From: Connor Holmes Date: Sat, 17 Aug 2024 12:49:36 -0400 Subject: [PATCH 1/3] Added an operation to convert an asymmetric polymatrix to a symmetric one. Added this functionality to initialization from sparse matrix (optional arg). --- _test/test_poly_matrix.py | 47 ++++++++++++++++++++++++++++++++++++++ poly_matrix/poly_matrix.py | 41 +++++++++++++++++++++++++++++---- 2 files changed, 84 insertions(+), 4 deletions(-) diff --git a/_test/test_poly_matrix.py b/_test/test_poly_matrix.py index bd84c07..1a79097 100644 --- a/_test/test_poly_matrix.py +++ b/_test/test_poly_matrix.py @@ -1,5 +1,7 @@ import numpy as np import pytest +import scipy.sparse as sp +import scipy.linalg as la from poly_matrix import PolyMatrix, sorted_dict @@ -448,6 +450,49 @@ def test_multiply(): np.testing.assert_allclose(S.get_matrix_dense(["x1", "x2"]), S_test) +def test_init_from_sparse(): + """Define polymatrix from sparse matrix. Test symmetrization.""" + # Get sparse matrix (already symmetric) + Q1 = get_Q() + var_dict = Q1.generate_variable_dict() + assert Q1.symmetric, ValueError("Starting matrix not symmetric. Test broken") + # Convert to sparse matrix + Q1_sparse = Q1.get_matrix() + # Convert back to polymatrix + Q2, _ = PolyMatrix.init_from_sparse(Q1_sparse, var_dict) + Q2_sparse = Q2.get_matrix() + # Test equality + np.testing.assert_allclose( + Q1_sparse.todense(), Q2_sparse.todense(), err_msg="Matrices not equal" + ) + + # Make matrix assymmetric by removing one block + Q2_sparse[4, 0] = 0.0 + Q2_sparse[4, 1] = 0.0 + Q2_sparse.eliminate_zeros() + Q2_mat = Q2_sparse.todense() + # Define new polymatrices + Q3_asym, _ = PolyMatrix.init_from_sparse(Q2_sparse, var_dict) + Q3_sym, _ = PolyMatrix.init_from_sparse(Q2_sparse, var_dict, symmetric=True) + # Check symmetry + Q3_asym_mat = Q3_asym.get_matrix(var_dict).todense() + assert not la.issymmetric(Q3_asym_mat), ValueError("Matrix should be assymmetric") + Q3_sym_mat = Q3_sym.get_matrix(var_dict).todense() + assert la.issymmetric(Q3_sym_mat), ValueError("Matrix should be symmetric") + # Check asymmetric matrix mapped properly + np.testing.assert_allclose( + Q2_mat, + Q3_asym_mat, + err_msg="Asymmetric matrix definition fail.", + ) + # Check that symmetry and get_matrix operations commute + np.testing.assert_allclose( + (Q2_mat + Q2_mat.T) / 2, + Q3_sym_mat, + err_msg="symmetry and matrix retrieval ops don't commute.", + ) + + if __name__ == "__main__": test_multiply() @@ -465,4 +510,6 @@ def test_multiply(): test_operations_simple() test_operations_advanced() + test_init_from_sparse() + print("all tests passed") diff --git a/poly_matrix/poly_matrix.py b/poly_matrix/poly_matrix.py index f641e7c..52e814d 100644 --- a/poly_matrix/poly_matrix.py +++ b/poly_matrix/poly_matrix.py @@ -102,8 +102,9 @@ def __init__(self, symmetric=True): self.symmetric = symmetric - # dictionary of form {variable-key : variable-size} + # dictionary of form {variable-key : variable-size} for enforcing # TODO(FD) consider replacing with NamedTuple + # consistent variable sizes for matrix blocks. self.variable_dict_i = {} self.variable_dict_j = {} @@ -126,7 +127,7 @@ def init_from_row_list(row_list, row_labels=None): return poly_vstack @staticmethod - def init_from_sparse(A, var_dict, unfold=False): + def init_from_sparse(A, var_dict, unfold=False, symmetric=False): """Construct polymatrix from sparse matrix (e.g. from learning method)""" self = PolyMatrix(symmetric=False) var_dict_augmented = augment(var_dict) @@ -150,10 +151,42 @@ def init_from_sparse(A, var_dict, unfold=False): if unfold: new_var_dict = {val[2]: 1 for val in var_dict_augmented.values()} return self, new_var_dict + + # Symmetrize if required + if symmetric: + self.make_symmetric() + return self, var_dict - def drop(self, variables): - for v in variables: + def make_symmetric(self): + """Convert a polymatrix from assymmetric to symmetric. + Performed in place.""" + if self.symmetric: + return + # Search through keys and identify + for key1 in self.matrix.keys(): + for key2 in self.matrix[key1].keys(): + # check if flipped order element exists (always holds for diagonals) + if key2 in self.matrix.keys() and key1 in self.matrix[key2].keys(): + # Symmetrize stored element + mat1 = self.matrix[key1][key2] + mat2 = self.matrix[key2][key1] + mat = (mat1 + mat2.T) / 2 + self.matrix[key1][key2] = mat + if not key1 == key2: + self.matrix[key2][key1] = mat.T + else: # If other side not populated, assume its zero + mat = self.matrix[key1][key2] / 2 + self.matrix[key1][key2] = mat + self.__setitem__((key2, key1), val=mat.T, symmetric=False) + # Align variable list and adjacency + self.variable_dict_j = self.variable_dict_i + self.adjacency_j = self.adjacency_i + # Set flag + self.symmetric = True + + def drop(self, variables_i): + for v in variables_i: if v in self.matrix: self.matrix.pop(v) if v in self.variable_dict_i: From 9c7b5623fb1e4623b6ca15ffa8c8eed8f8f25ec8 Mon Sep 17 00:00:00 2001 From: Connor Holmes Date: Sat, 17 Aug 2024 17:38:56 -0400 Subject: [PATCH 2/3] updated expression generation so that homogeneous variable can be excluded if desired. --- poly_matrix/poly_matrix.py | 32 ++++++++++++++++++++------------ 1 file changed, 20 insertions(+), 12 deletions(-) diff --git a/poly_matrix/poly_matrix.py b/poly_matrix/poly_matrix.py index 52e814d..3b2e116 100644 --- a/poly_matrix/poly_matrix.py +++ b/poly_matrix/poly_matrix.py @@ -673,7 +673,7 @@ def get_block_matrices(self, key_list=None): blocks.append(np.zeros((i_size, j_size))) return blocks - def get_expr(self, variables=None): + def get_expr(self, variables=None, homog=None): if not self.symmetric: warnings.warn( @@ -708,19 +708,27 @@ def get_expr(self, variables=None): if i > 0: sub_expr += " + " if diag: - sub_expr += ( - "{:.3f}".format(vals[i]) - + "*" - + ":".join([key1, str(row[i])]) - + "^2" - ) + if key1 == homog: + sub_expr += "{:.3f}".format(vals[i]) + else: + sub_expr += ( + "{:.3f}".format(vals[i]) + + "*" + + ":".join([key1, str(row[i])]) + + "^2" + ) else: + if key1 == homog: + key1str = "" + else: + key1str = "*" + ":".join([key1, str(row[i])]) + if key2 == homog: + key2str = "" + else: + key2str = "*" + ":".join([key2, str(col[i])]) + sub_expr += ( - "{:.3f}".format(vals[i] * 2) - + "*" - + ":".join([key1, str(row[i])]) - + "*" - + ":".join([key2, str(col[i])]) + "{:.3f}".format(vals[i] * 2) + key1str + key2str ) if first_expr: expr += sub_expr + "\n" From 722d2c64cf8173d838035829dcbf2a7aa8db3610 Mon Sep 17 00:00:00 2001 From: Connor Holmes Date: Tue, 24 Sep 2024 11:27:53 -0400 Subject: [PATCH 3/3] Remove obsolete adjacency definition in make_symmetric method --- poly_matrix/poly_matrix.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/poly_matrix/poly_matrix.py b/poly_matrix/poly_matrix.py index 3b2e116..ec52f2b 100644 --- a/poly_matrix/poly_matrix.py +++ b/poly_matrix/poly_matrix.py @@ -179,9 +179,8 @@ def make_symmetric(self): mat = self.matrix[key1][key2] / 2 self.matrix[key1][key2] = mat self.__setitem__((key2, key1), val=mat.T, symmetric=False) - # Align variable list and adjacency + # Align variable list self.variable_dict_j = self.variable_dict_i - self.adjacency_j = self.adjacency_i # Set flag self.symmetric = True