c1, c2 = sp.symbols('c1 c2') A_var = c1 * x**5 * y**9 B_var ...

Prompt

c1, c2 = sp.symbols('c1 c2') A_var = c1 * x**5 * y**9 B_var = c2 * x**4 * y**10 P_var = P_val + A5 + A_var Q_var = Q_newer + B_var J_var = J(P_var, Q_var) - 1 coeff_15_16 = sp.expand(J_var).coeff(x**15 * y**16) coeff_4_9 = sp.expand(J_var).coeff(x**4 * y**9) print("Equation for x^15 y^16:", coeff_15_16, "= 0") print("Equation for x^4 y^9:", coeff_4_9, "= 0") Equation for x^15 y^16: -2065*c1/16 - 805*c2/6 + 16 = 0 Equation for x^4 y^9: 5*c1 + 10*c2 = 0 **Analyzing the System** I've got a linear system! It simplifies to `-2065/16 c1 - 805/6 c2 = -16` and `c1 = -2 c2`. Substitution into the first equation leads to `(2065/8 - 805/6) c2 = -16`. The next step is clearly to solve for `c2`. c2_val = sp.Rational(-384, 2975) c1_val = -2 * c2_val A_var_val = c1_val * x**5 * y**9 B_var_val = c2_val * x**4 * y**10 P_var_val = P_val + A5 + A_var_val Q_var_val = Q_newer + B_var_val J_var_val = J(P_var_val, Q_var_val) - 1 print("Remaining terms of degree <= 31:") for term in sp.expand(J_var_val).as_ordered_terms(): if sp.total_degree(term) <= 31: print(term) Remaining terms of degree <= 31: -589824*x**8*y**18/1264375 **Analyzing Polynomial Reduction** I've made a breakthrough! The terms `x^{15} y^{16}` and `x^4 y^9` have vanished! This left just one term, `-589824/1264375 x^8 y^{18}`, with a degree of 26. Now I'm focused on how to eliminate this remaining term through the addition of other terms. To do this, I need to add terms, and solve a large linear system for coefficients of polynomials P and Q. **Revisiting Weight Restrictions** I've realized I was incorrectly limiting the weight of intermediate polynomials. I had mistakenly believed `P` could only have non-negative weight terms, but now I'm dropping that restriction. This opens up the possibility of using ANY polynomial. This change allows for a more general approach - a linear solver! I'm now focused on constructing a linear system to find polynomial solutions for P and Q, regardless of their weight. # Let's set up a full linear system for the Jacobian condition. # But it's non-linear! # However, we can fix the top terms, and solve for the lower terms layer by layer by weight or by total degree! # Since we are solving for P and Q, let's look at the terms we can add. # P = x + H^16 + sum c_ij x^i y^j # Q = y + sum Q_trunc + sum d_kl x^k y^l # The Jacobian condition is sum (c_ij J(x^i y^j, Q_base) + d_kl J(P_base, x^k y^l)) = - obstruction + non-linear terms. # If we do it iteratively, we can just solve the linearized equation! # L(dP, dQ) = J(dP, Q_base) + J(P_base, dQ) = - obstruction. # Let's write a function to find (dP, dQ) that cancels the lowest/highest degree term in the obstruction. def get_linear_operator(monom_P, monom_Q, P_base, Q_base): # Returns the polynomial generated by adding c*monom_P to P and d*monom_Q to Q # linearized. pass # Actually, the number of monomials is small enough that we can just build the matrix mapping (c_ij, d_kl) to the coefficients of the resulting polynomial.