allowing symbolic powers

qbee/polynomialization.py

@@ -230,10 +230,11 @@ def _try_convert_to_polynomial(self, expr: sp.Expr) -> sp.Expr | None:
                 replaced_non_pow_lazy = expr.subs(
                     {v: k for k, v in self._substitution_equations.items() if not v.is_Pow})
             replaced_pow_lazy = replaced_non_pow_lazy \
+                .replace(sp.Pow, self._replace_symbolic_pow) \
                 .replace(sp.Pow, self._replace_irrational_pow) \
                 .replace(sp.Pow, self._replace_negative_integer_pow)
         except Exception as e:
-            warnings.warn("Substituting new variables failed to produce the expected calculations,"
+            earnings.warn("Substituting new variables failed to produce the expected calculations,"
                           " so the calculation was done in an alternative mode. "
                           "If you see this message, please let us know in the Issues: "
                           "https://github.com/AndreyBychkov/QBee/issues", RuntimeWarning)
@@ -246,6 +247,7 @@ def _try_convert_to_polynomial(self, expr: sp.Expr) -> sp.Expr | None:
         # eager evaluation
         replaced_non_pow = expr.subs({v: k for k, v in self._substitution_equations.items() if not v.is_Pow})
         replaced_pow = replaced_non_pow \
+            .replace(sp.Pow, self._replace_symbolic_pow) \
             .replace(sp.Pow, self._replace_irrational_pow) \
             .replace(sp.Pow, self._replace_negative_integer_pow)
         return replaced_pow if self._is_expr_polynomial(replaced_pow) else None
@@ -263,6 +265,14 @@ def _replace_irrational_pow(self, base, exp):
             return new_vars[0] if len(new_vars) != 0 else base ** exp
         return base ** exp
 
+    def _replace_symbolic_pow(self, base, exp):
+        if not (exp.is_Float or exp.is_Integer):
+            for k, v in self._substitution_equations.items():
+                if v.is_Pow and v.base == base:
+                    if sp.simplify(exp / v.exp).is_Integer:
+                        return k ** sp.simplify(exp / v.exp)
+        return base ** exp
+
     def _is_expr_polynomial(self, expr: sp.Expr):
         if not self.variables.laurent:  # performance optimization
             return expr.is_polynomial(*self.variables.state, *self.variables.input)
@@ -457,8 +467,17 @@ def available_substitutions(system: EquationSystem) -> set[sp.Expr]:
     return filter_laurent_monoms(system, unused_subs)
 
 
-def find_pow(values: Collection[sp.Expr], pow: sp.Pow) -> sp.Pow | None:
-    all_matches = [v for v in values if v.is_Pow and v.base == pow.base and math.isclose(v.exp, pow.exp)]
+def find_pow(values: Collection[sp.Expr], pw: sp.Pow) -> sp.Pow | None:
+    all_matches = []
+    numtypes = (sp.core.numbers.Float, sp.core.numbers.Integer, float, int)
+    for v in values:
+        if v.is_Pow and v.base == pw.base:
+            if type(v.exp) in numtypes and type(pw.exp) in numtypes:
+                if math.isclose(v.exp, pw.exp):
+                    all_matches.append(v)
+            else:
+                if sp.simplify(v.exp - pw.exp) == 0:
+                    all_matches.append(v)
     if all_matches:
         return all_matches[0]
     return None

tests/test_polynomialization.py

@@ -94,6 +94,30 @@ def test_nested_functions():
     res = polynomialize(system)
     assert len(res) == 4
 
+def test_parametric_pow():
+    x = functions("x")
+    a = parameters("a")
+    system = [(x, x**a)]
+    res = polynomialize(system)
+    assert len(res) == 2
+
+    x = functions("x", laurent=False)
+    res = polynomialize([(x, x**a)])
+    assert len(res) == 3
+
+    res = polynomialize([(x, x**(-a))])
+    assert len(res) == 3
+
+def test_hill():
+    x = functions("x")
+    a, b = parameters("a b")
+    system = [(x, 1 / (a + x**b))]
+    res = polynomialize(system)
+    assert len(res) == 3
+
+    x = functions("x", laurent=False)
+    res = polynomialize([(x, 1 / (a + x**b))])
+    assert len(res) == 4
 
 def test_parameter():
     x, y = functions("x, y", laurent=False)