add comments and remove unused code in denserepr module

2 files changed, 41 insertions, 68 deletions

diff --git a/polymatrix/denserepr/from_.py b/polymatrix/denserepr/from_.py
index aa2e015..7936c97 100644
--- a/polymatrix/denserepr/from_.py
+++ b/polymatrix/denserepr/from_.py
@@ -23,6 +23,25 @@ def from_polymatrix(
 ) -> StateMonadMixin[
     ExpressionState, tuple[tuple[tuple[np.ndarray, ...], ...], tuple[int, ...]]
 ]:
+    """
+    Converts a tuple of polynomial vectors `expressions` into a (semi-dense) matrix representation.
+
+    Given single polynomial vector
+
+        expr = [[x**2 + y], [1 + x*y + y**2]],
+
+    its matrix representation is given by three matrices
+
+        1. np.array type matrix [[0], [1]] containing the constant part
+        2. np.array type matrix [[0, 1], [0, 0]] containing the linear part 
+           w.r.t the monomial vector Z=[x, y]
+        3. scipy.sparse type matrix [[1, 0, 0, 0], [0, 0.5, 0.5, 1]] containing the quadratic part
+           w.r.t the monomial vector Z=[x**2, x*y, x*y, y**2].
+
+    Instead of providing only a single polynomial vector, a tuple of polynomial vectors can be provided.
+    This approach ensures consistency in the vector of monomial accross all matrix representations.
+    """
+
     if isinstance(expressions, Expression):
         expressions = (expressions,)
 
@@ -78,6 +97,10 @@ def from_polymatrix(
                     f"{duplicates=}. Make sure you give a unique name for each variables."
                 )
 
+        # This dictionary maps the expression variable index (saved in the expression state) 
+        # to the variable index of the dense representation.
+        # Assume state.indices = {'x': (2,3), 'y': (1,2), 'z': (0,1)} and variable = [[2], [1]],
+        # then variable_index_map = {2: 0, 1: 1}.
         variable_index_map = {old: new for new, old in enumerate(sorted_variable_index)}
 
         n_param = len(sorted_variable_index)
@@ -92,6 +115,7 @@ def from_polymatrix(
                     n_param=n_param,
                 )
 
+                # MS: ignores any columns, better would be to check that there are no columns
                 for row in range(n_row):
                     polymatrix_terms = polymatrix.get_poly(row, 0)
 
@@ -117,12 +141,21 @@ def from_polymatrix(
                                     for _ in range(count):
                                         yield index
 
+                            # converts the monomial x**2 = ((2, 2),) to (2, 2)
+                            # converts the monomial x*y = ((2, 1), (1, 1)) to (2, 1)
+                            # converts the monomial y**2 = ((1, 2),) to (1, 1)
                             new_variable_indices = tuple(gen_new_monomial())
 
+                            # converts (2, 2) to (0,)
+                            # converts (2, 1) to (1, 2)
+                            # converts (1, 2) to (3,)
+                            # the cols correspond to the column of the dense matrix w.r.t. Z
                             cols = variable_indices_to_column_index(
                                 n_param, new_variable_indices
                             )
 
+                            # the monomial x*y is mapped to two columns, hence we divide its value by 2:
+                            # x*y = [0, 0.5, 0.5, 0] @ [x**2, x*y, x*y, y**2].T
                             col_value = value / len(cols)
 
                             for col in cols:
@@ -133,47 +166,8 @@ def from_polymatrix(
 
         underlying_matrices = tuple(gen_numpy_matrices())
 
-        def gen_auxillary_equations():
-            for key, monomial_terms in state.auxillary_equations.items():
-                if key in sorted_variable_index:
-                    yield key, monomial_terms
-
-        auxillary_equations = tuple(gen_auxillary_equations())
-
-        n_row = len(auxillary_equations)
-
-        if n_row == 0:
-            auxillary_matrix_equations = None
-
-        else:
-            buffer = DenseReprBufferImpl(
-                data={},
-                n_row=n_row,
-                n_param=n_param,
-            )
-
-            for row, (key, monomial_terms) in enumerate(auxillary_equations):
-                for monomial, value in monomial_terms.items():
-                    new_monomial = tuple(
-                        variable_index_map[var]
-                        for var, count in monomial
-                        for _ in range(count)
-                    )
-
-                    cols = variable_indices_to_column_index(n_param, new_monomial)
-
-                    col_value = value / len(cols)
-
-                    for col in cols:
-                        buffer.add(
-                            row, col, sum(count for _, count in monomial), col_value
-                        )
-
-            auxillary_matrix_equations = buffer
-
         result = DenseReprImpl(
             data=underlying_matrices,
-            aux_data=auxillary_matrix_equations,
             variable_mapping=sorted_variable_index,
             state=state,
         )
diff --git a/polymatrix/denserepr/impl.py b/polymatrix/denserepr/impl.py
index 294a396..9f53c8b 100644
--- a/polymatrix/denserepr/impl.py
+++ b/polymatrix/denserepr/impl.py
@@ -46,45 +46,24 @@ class DenseReprBufferImpl:
 # NP: as discussed in meeting #6, this is not actually dense
 # NP: also representation here does not mean much, this is actually a series of
 # NP: matrices from a cone, consider changing name to reflet this fact.
-
 # NP: also this really needs documentation on how it is indexed, it is _not_
 # NP: intuitive from the method names
 @dataclasses.dataclass
 class DenseReprImpl:
     data: tuple[DenseReprBufferImpl, ...]
-    aux_data: typing.Optional[DenseReprBufferImpl]
     variable_mapping: tuple[int, ...]
     state: ExpressionState
 
-    # NP: as discussed in meeting #5 the "equations" are not equations
-    # NP: rename to something more meaningful
-    def merge_matrix_equations(self):
-        # FIXME: this function is never used and shadowed by the other gen_matrices below
-        def gen_matrices(index: int):
-            for equations in self.data:
-                if index < len(equations):
-                    yield equations[index]
-
-            if index < len(self.aux_data):
-                yield self.aux_data[index]
-
-        indices = set(
-            key
-            for equations in self.data + (self.aux_data,)
-            for key in equations.keys()
-        )
-
-        def gen_matrices():
-            for index in indices:
-                if index <= 1:
-                    yield index, np.vstack(tuple(gen_matrices(index)))
-                else:
-                    yield index, scipy.sparse.vstack(tuple(gen_matrices(index)))
-
-        return dict(gen_matrices())
-
     # NP: get variable value out of solution vector
     def get_value(self, variable, value): # NP: variable index (state), value result from optimiz
+        """
+        Given 
+            variable = 'x'
+            value = np.array((0.1, 0.2, 0.3, 0.4, 0.5))
+            state = {'x': (1, 3)}
+        the function returns np.array((0.2, 0.3)).
+        """
+        
         variable_indices = get_variable_indices_from_variable(self.state, variable)[1]
 
         def gen_value_index():