diff options
| -rw-r--r-- | sumofsquares/canon.py | 15 | ||||
| -rw-r--r-- | sumofsquares/solver/cvxopt.py | 1 |
2 files changed, 8 insertions, 8 deletions
diff --git a/sumofsquares/canon.py b/sumofsquares/canon.py index 02c9266..8c4b541 100644 --- a/sumofsquares/canon.py +++ b/sumofsquares/canon.py @@ -114,7 +114,7 @@ class PutinarPSatz(Canonicalization): new_constraints: list[Constraint] = list(constraints) for i, constr in enumerate(to_process): new_constr = constr.expression - for domain_poly in constr.domain.polynomials: + for j, domain_poly in enumerate(constr.domain.polynomials): # Make a multiplier polynomial for a domain polynomial of a constraint. # # Since Putinar's PSatz does not explicity tell us the order of the @@ -140,20 +140,19 @@ class PutinarPSatz(Canonicalization): if d.read(state).scalar().constant() % 2 != 0: d += 1 - x = poly.v_stack((1,) + tuple( - init_variable_expr(v, state.get_shape(v)) - for v in prob.polynomial_variables)) + x = poly.v_stack((1,) + tuple(init_variable_expr(v, state.get_shape(v)) + for v in prob.polynomial_variables)) - # FIXME: need to check that there is not a \gamma_{i} already! + # FIXME: need to check that there is not a \gamma_i,j already! # TODO: deterministically generate unique names # TODO: better documentation # Generate monomials of up to degree d and take only monomials # that are in half of the newton polytope - monomials = x.combinations(poly.arange(d + 1)).half_newton_polytope(x) - c = opt_variable_from_name(rf"\gamma_{i}", shape=monomials.shape) + monomials = x.combinations(poly.arange(d + 1)) # .half_newton_polytope(x) + gamma = opt_variable_from_name(rf"\gamma_{i},{j}", shape=monomials.shape) - multiplier = c.T @ monomials + multiplier = gamma.T @ monomials # Subtract multiplier from expression and impose that it is also SOS new_constr -= multiplier * domain_poly diff --git a/sumofsquares/solver/cvxopt.py b/sumofsquares/solver/cvxopt.py index 0692787..cf0dc4e 100644 --- a/sumofsquares/solver/cvxopt.py +++ b/sumofsquares/solver/cvxopt.py @@ -135,6 +135,7 @@ def solve_cone(prob: ConicProblem, verbose: bool = False, info = cvxopt.solvers.coneqp(P=P, q=q, G=G, h=h, A=A, b=b, dims=dims, *args, **kwargs) else: + assert P is None, "Solving QP with LP solver! Something has gone very wrong" info = cvxopt.solvers.conelp(c=q, G=G, h=h, A=A, b=b, *args, **kwargs) |