coniclifts is sageopt’s backend

Sageopt does not solve SAGE relaxations on its own; it relies on third party convex optimization solvers, such as ECOS or MOSEK. These solvers require input in very specific standard-forms. Coniclifts provides abstractions that allow us to state SAGE relaxations in high-level syntax, and manage interactions with these low-level solvers.

Overview

Coniclifts employs abstractions that are similar to those in CVXPY. For example, it is possible to construct and solve a linear program using coniclifts:

import sageopt.coniclifts as cl
import numpy as np
# random problem data
G = np.random.randn(3, 6)
h = A @ np.random.rand(6)
c = np.random.rand(6)
# input to coniclift's Problem constructor
x = cl.Variable(shape=(6,))
constrs = [0 <= x, G @ x == h]
objective_expression = c @ x
prob = cl.Problem(cl.MIN, objective_expression, constrs)
prob.solve(solver='ECOS', verbose=False)
x_opt = x.value

In the example above, constrs consists of two coniclifts Constraint objects. Both of these objects were constructed with operator overloading (<= and ==), and the second Constraint made use of operator-overloaded matrix multiplication (@). If you check the datatype of the object y = A @ x, you would find that y is a coniclifts Expression. Expression objects track functions of Variable objects.

Coniclifts has support for nonlinear convex constraints. The most important of these constraints are specified as “set membership”, rather than elementwise inequalities. For example, we have a PrimalSageCone class, which is used to construct some of the SAGE relaxations in the sageopt.relaxations package. Here is a concrete demonstration

alpha = np.array([[0, 0],
                  [1, 0],
                  [0, 1],
                  [1, 1],
                  [0.5, 0],
                  [0, 0.5]])
c = np.array([0, 3, 2, 1, -4, -2])
# alpha, c define a signomial in the usual way
gamma = cl.Variable(shape=(), name='gamma')
c_expr = c.copy()
c_expr[0] -= gamma   # shift the constant term by -gamma
constr = cl.PrimalSageCone(c_expr, alpha, X=None, name='example_constraint')
prob = Problem(cl.MAX, gamma, [constr])
# find largest gamma so shifted signomial is nonnegative
status, val = prob.solve()

By solving the problem described above, we have that val is a lower bound on the signomial which takes values lambda x: c @ np.exp(alpha @ x). There are only a few set-membership constraints currently implemented in coniclifts; they can all be found in sageopt.coniclifts.constraints.set_membership.

Coniclifts contains operators which represent functions applied to Expression objects. The module sageopt/coniclifts/operators/affine.py contains all affine array-manipulation operators you would use on a numpy ndarray [link]. This includes linear algebra operators such as kron or trace, reshaping operators such as hstack or tile, and array-creation routines like diag or triu. These functions behave in identical ways to their numpy counterparts, because Expression objects are actually a custom subclass of numpy’s ndarray datatype.

Coniclifts also has a small selection of nonlinear operators: weighted_sum_exp, vector2norm, and relent. It is easy to add more nonlinear operators, but these three suffice for the internal uses currently found in sageopt. Here is a concrete example

alpha = np.array([[1, 0],
                  [0, 1],
                  [1, 1],
                  [0.5, 0],
                  [0, 0.5]])
c = np.array([3, 2, 1, 4, 2])
x = cl.Variable(shape=(2,), name='x')
cons = [cl.weighted_sum_exp(c, alpha @ x) <= 1]
obj = -x[0] - 2*x[1]
prob = Problem(cl.MIN, obj, cons)
prob.solve()
x_expect = np.array([-4.93083, -2.73838])
x_actual = x.value
print(np.allclose(x_expect, x_actual, atol=1e-4))

If you run the code above, you should find that it prints True.

It is important to note that nonlinear operators are not allowed in the objective function. So if you want to minimize a nonlinear convex function given by expr, you need to create an auxiliary variable such as t = Variable(shape=(1,)), add the constraint expr <= t, and set the objective to minimize t. Here is an example of a constrained least-squares problem we solve for solution recovery in dual SAGE relaxations

A, b, K = con.X.A, con.X.b, con.X.K  # con is a DualSageCone instance
log_v = np.log(con.v.value)
A = np.asarray(A)
x = cl.Variable(shape=(A.shape[1],))
t = cl.Variable(shape=(1,))
cons = [cl.vector2norm(log_v - alpha @ x[:con.n]) <= t,
        cl.PrimalProductCone(A @ x + b, K)]
prob = cl.Problem(cl.MIN, t, cons)

The example above also alludes to a useful set-membership constraint, called PrimalProductCone. Follow this link for source code of set-membership constraints available in coniclifts.

The Variable class

class sageopt.coniclifts.Variable

An abstraction for a symbol appearing in constraint sets, or optimization problems.

Variable objects are a custom subclass of numpy ndarrays.

Parameters
  • shape (tuple) – The dimensions of the Variable object. Defaults to shape=().

  • name (str) – A string which should uniquely identify this Variable object in all models where it appears. Ideally, this string should be human-readable. Defaults to 'unnamed_var_N', where N is an integer.

  • var_properties (list of str) – Currently, the only accepted forms of this argument are the empty list (in which case an unstructured Variable is returned), or a list containing the string 'symmetric' (in which case a symmetric matrix Variable is returned).

Examples

The symbol you use in the Python interpreter does not need to match the “name” of a Variable.

x = Variable(shape=(3,), name='my_name')

A Variable object can take on any dimension that a numpy ndarray could take on.

y = Variable(shape=(10,4,1,2), name='strange_shaped_var')

Notes

Upon construction, Variable objects are “proper”. If you index into them, they are still considered Variable objects, but they no longer contain information about all of their components. A Variable object’s name field only uniquely determines the “proper” version of that Variable. If v.is_proper() == False, then it should be possible to recover the original Variable object with original_v = v.base.

x = Variable(shape=(3,), name='x')
print(type(x))  # sageopt.coniclifts.base.Variable
print(x.is_proper())  # True

y = x[1:]
print(type(y))  # sageopt.coniclifts.base.Variable
print(y.is_proper())  # False
print(x.name == y.name)  # True
print(id(x) == id(y.base))  # True; these exist at the same place in memory.
property generation

An internally-maintained index. Variable objects of different “generation” cannot participate in a common optimization problem.

property name

A string which should uniquely identify this object in all models where it appears, provided self.is_proper() == True.

property scalar_variable_ids

Each component of this Variable object (i.e. each “scalar variable”) contains an index which uniquely identifies it in all models where this Variable appears. Return the list of these indices.

scalar_variables()

Return a list of all ScalarVariables appearing in this Expression.

property value

An ndarray containing numeric entries, of shape equal to self.shape. This is the result of the most recent call to set_scalar_variables.

The Problem class

Detailed documentation for the Problem class is given below. Most users will only interact with a few aspects of Problem objects. On a first read, it should be enough just to skim the documentation for this class. If you want to understand all the attributes of a Problem object, you will need to read The Expression system and The compiler interface.

class sageopt.coniclifts.Problem(objective_sense, objective_expr, constraints)

A representation for a convex optimization problem. When this Problem object is constructed, the constraints are immediately compiled into a flattened conic representation. This compilation process may take some time for large problems.

Parameters
  • objective_sense (str) – Either coniclifts.MIN or coniclifts.MAX

  • objective_expr (Expression) – The function to minimize or maximize. ScalarExpressions and real numeric types are also accepted. The final (cast) expression must be linear in its Variables.

  • constraints (list of coniclifts.Constraint) – Variables appearing in objective_expr are subject to Constraint objects in this list.

objective_sense

The value passed by the user to the constructor.

Type

str

objective_expr

The value passed by the user to the constructor.

Type

Expression

constraints

The value passed by the user to the constructor.

Type

list of coniclifts.Constraint

A

The matrix in the flattened conic representation of the feasible set.

Type

CSC-format sparse matrix

b

The offset vector in the flattened conic representation of the feasible set.

Type

ndarray

K

The cartesian product of these cones (in order) defines the convex cone appearing in the flattened conic representation of the feasible set.

Type

list of coniclifts.Cone

all_variables

All Variable objects needed to represent the feasible set. This includes user-declared Variables, and Variables which were required to write the problem in terms of coniclifts primitives. It is recommended to access this list by calling self.variables().

Type

list of coniclifts.Variable

variable_map

A map from a Variable’s name field to a numpy array. If myvar is a coniclifts Variable appearing in the system defined by constraints, then a point x satisfying \(A x + b \in K\) maps to a feasible value for myvar by

x0 = np.hstack([x, 0])
myvar_val = x0[variable_map[myvar.name]]

In particular, we guarantee myvar.shape == variable_map[myvar.name].shape. Augmenting x by zero to create x0 reflects a convention that if a component of a Variable does not affect the constraints, that component is automatically assigned the value zero.

Type

Dict[str,ndarray]

variable_values

A map from a Variable’s name field to a numpy array, containing a feasible value for that Variable for this problem’s constraints.

Type

Dict[str, ndarray]

solver_apply_data

A place to store metadata during a call to self.solve(cache_apply_data=True).

Type

Dict[str,dict]

solver_raw_output

A place to store metadata during a call to self.solve(cache_raw_output=True).

Type

Dict[str,dict]

metadata

A place for users to store metadata produced when constructing this Problem.

Type

dict

value
The objective value from the last call to self.solve. Can be a float,

np.inf, -np.inf, or np.NaN.

Type

float

status

The problem status from the last call to self.solve. Either coniclifts.SOLVED or coniclifts.INACCURATE or coniclifts.FAILED.

Type

str

Notes

Problem status being coniclifts.SOLVED does not mean that the decision variables have been assigned specific values. It only means that the solver returned a normal status code (i.e. that the solver didn’t run into numerical difficulties). If a solver return code indicates the problem is “infeasible” or “unbounded”, we still set self.status = coniclifts.SOLVED. Therefore it is important to check that self.status == coniclifts.SOLVED and that -np.inf < self.value < np.inf before accessing a Variable’s value.

solve(solver=None, **kwargs)
Parameters

solver (str) – None or 'MOSEK' or 'ECOS'. Defaults to the first installed solver found in Problem._SOLVER_ORDER_.

Returns

  • status (str) – Either coniclifts.SOLVED or coniclifts.INACCURATE or coniclifts.FAILED. Refer to the Notes of the Problem class for the meanings of these values.

  • value (float) – The optimal objective value reported by the solver. Can be a float, -np.inf, np.inf, or np.NaN.

Notes

Keyword arguments.

cache_apply_data

variables()

Return a shallow copy of self.all_variables.

This function is provided to match the syntax of CVXPY Problem objects.

The Expression system

Coniclifts is built around a few core ideas, including …

  • transparency in the compilation process,

  • ease-of-extension for experts in convex optimization,

  • no dependence on a C or C++ backend,

  • full compatibility with numpy.

In order to achieve full compatibility with numpy, coniclifts takes an elementwise approach to symbolic expressions. Specifically, coniclifts begins with a few simple abstractions for scalar-valued symbolic expressions, and wraps those abstractions in a custom subclass of numpy’s ndarray. The coniclifts abstractions for scalar-valued symbolic expressions are as follows:

  • A ScalarExpression class represents scalar-valued affine functions of certain irreducible primatives. ScalarExpressions are operator-overloaded to support +, -, and *. This allows ndarrays of ScalarExpressions to fall back on many functions which are implemented for numeric ndarrays.

  • An abstract ScalarAtom class specifies the behavior of the irreducible primitives in ScalarExpressions. The ScalarAtom class immediately specializes into ScalarVariables (far and away the most important ScalarAtom), and another abstract class, called NonlinearScalarAtom. NonlinearScalarAtoms are implemented on a case-by-case basis, but include such things as the exponential function and the vector 2-norm.

We ask interested users to refer to the source code for additional information on ScalarExpressions and ScalarAtoms. For most people, all you need to work with is the Expression class.

class sageopt.coniclifts.Expression

An Expression is an ndarray whose entries are ScalarExpressions. Variable objects are a special case of Expression objects.

Construction

Expression objects can be constructed from ndarrays of numeric types, or ndarrays containing ScalarExpressions. In both cases, the construction syntax is

expr = Expression(existing_array)

Arithmetic operator overloading

Operations +, -, / and * work in the exact same way as for numpy arrays.

Expressions do not allow exponentiation (**).

Expressions overload @ (a.k.a. __matmul__) in a way that is consistent with numpy, but only for arguments which have up to two dimensions.

Constraint-based operator overloading

Operations <= and >= produce elementwise inequality constraints.

The operation == produces an elementwise equality constraint.

The operations >> and << produce linear matrix inequality constraints.

static are_equivalent(expr1, expr2, rtol=1e-05, atol=1e-08)

Perform a check that expr1 and expr2 are symbolically equivalent, in the sense of affine operators applied to ScalarAtoms. The equivalence is up to numerical tolerance in the sense of np.allclose.

Parameters
  • expr1 (Expression) –

  • expr2 (Expression) –

  • rtol (float) – relative numerical tolerance

  • atol (float) – absolute numerical tolerance

Returns

Return type

True if the Expressions can be verified as symbolically equivalent. False otherwise.

Notes

The behavior of this function is conservative. If self contains a mix of ScalarAtoms (e.g. ScalarVariables and NonlinearScalarAtoms), then this function might return False even when expr1 and expr2 are equivalent. This is due to nondeterministic behavior of Expression.factor in such situations.

factor()

Returns a tuple (A, x, B).

A is a tensor of one order higher than the current Expression object, i.e. A.ndim == self.ndim + 1. The dimensions of A and self agree up until self.ndim, i.e. A.shape[:-1] == self.shape.

x is a list of ScalarAtom objects, with len(x) == A.shape[-1].

B is a numpy array of the same shape as self.

The purpose of this function is to enable faster matrix multiplications of Expression objects. The idea is that if you tensor-contract A along its final dimension according to x, and then add B, you recover this Expression.

is_affine()

True if the Expression is an affine function of Variables within its scope.

is_concave()

Return an ndarray of booleans. For a fixed component index, the value of the returned array indicates if that component of the current Expression is a concave function of Variables within its scope.

is_constant()

True if no Variables appear in this Expression.

is_convex()

Return an ndarray of booleans. For a fixed component index, the value of the returned array indicates if that component of the current Expression is a convex function of Variables within its scope.

scalar_atoms()

Return a list of all ScalarAtoms appearing in this Expression.

scalar_variables()

Return a list of all ScalarVariables appearing in this Expression.

property value

An ndarray containing numeric entries, of shape equal to self.shape. This is the result of propagating the value of ScalarVariable objects through the symbolic operations tracked by this Expression.

variables()

Return a list of all Variable objects appearing in this Expression. You can assume that all Variable objects will be “proper”.

SAGE constraint classes

Coniclifts provides direct implementations of the primal and dual signomial SAGE cones. The implementation details between ordinary-SAGE and conditional-SAGE versions of the primal and dual cones are abstracted away by the classes

These classes have virtually identical constructors and public attributes. In particular, both classes’ constructors require an argument X, which can be None or a SigDomain. Ordinary SAGE constraints are obtained by setting X=None. Conditional SAGE constraints assume the feasible set represented by induced by X is feasible, and it is the user’s responsibility to ensure this is the case. The main difference in these classes’ attributes is that PrimalSageCone instances have a dict called age_vectors (which represent the certificates of nonnegativity) and that DualSageCone instances have a dict called mu_vars (which are useful for solution recovery in SAGE relaxations of signomial programs).

The PrimalSageCone class performs a very efficient dimension-reduction procedure by analyzing the signs of the provided vector c. The details of the reduction are described in Corollary 5 of MCW2019. At present, coniclifts does not provide a means to track the signs of Variable objects, and so this reduction is limited to indices i where c[i] is constant. This feature can optionally be carried over to DualSageCone objects, if the user provides a keyword argument c to the DualSageCone constructor.

The PrimalSageCone and DualSageCone classes can perform a more extensive presolve phase to eliminate trivial AGE cones (those which reduce to the nonnegative orthant). By default, this presolve ability is turned off. This default can be changed by calling sageopt.coniclifts.presolve_trivial_age_cones(True). The computational cost of this presolve is borne when the constraint is constructed, and scales linearly in the dimension of the SAGE constraint (equal to constr.alpha.shape[0]). The cost of this presolve can be mitigated by recycling covers = constr.ech.expcovers from one call of a constraint constructor to the next.

class sageopt.coniclifts.PrimalSageCone(c, alpha, X, name, **kwargs)

Represent the constraint “\(c \in C_{\mathrm{SAGE}}(\alpha, X)\)”.

Parameters
  • c (Expression) – The vector subject to this PrimalSageCone constraint.

  • alpha (ndarray) – The rows of alpha are the exponents defining this primal SAGE cone.

  • X (SigDomain or None) – If None, then this constraint represents a primal ordinary-SAGE cone.

  • name (str) – Uniquely identifies this Constraint in the model where it appears. Serves as a suffix for the name of any auxiliary Variable created when compiling to the coniclifts standard.

Other Parameters

covers (Dict[int, ndarray]) – covers[i] indicates which indices j have alpha[j,:] participate in the i-th AGE cone. Automatically constructed in a presolve phase, if not provided.

alpha

The rows of alpha are the exponents defining this primal SAGE cone.

Type

ndarray

c

The vector subject to this PrimalSageCone constraint.

Type

Expression

age_vectors

If all Variable objects in the scope of this constraint are assigned feasible, values, then we should have age_vectors[i].value in the i-th AGE cone with respect to alpha, and c.value == sum([av.value for av in age_vectors.values()], axis=0).

Type

Dict[int, Expression]

X

If None, then this constraint represents a primal ordinary-SAGE cone.

Type

SigDomain or None

ech

A simple wrapper around the constructor argument covers. Manages validation of covers when provided, and manages construction of covers when a user does not provide it.

Type

ExpCoverHelper

name

Uniquely identifies this Constraint in the model where it appears. Serves as a suffix for the name of any auxiliary Variable created when compiling to the coniclifts standard.

Type

str

Notes

The constructor can raise a RuntimeError if the constraint is deemed infeasible.

class sageopt.coniclifts.DualSageCone(v, alpha, X, name, **kwargs)

Represent the constraint “\(v \in C_{\mathrm{SAGE}}(\alpha, X)^{\dagger}\)”.

Parameters
  • v (Expression) – The vector subject to the dual SAGE-cone constraint.

  • alpha (ndarray) – The matrix of exponent vectors defining the SAGE cone; alpha.shape[0] == v.size.

  • X (SigDomain or None) – If None, then this constraint represents a dual ordinary-SAGE cone.

  • name (str) – Uniquely identifies this Constraint in the model where it appears. Serves as a suffix for the name of any auxiliary Variable created when compiling to the coniclifts standard.

Other Parameters
  • covers (Dict[int, ndarray]) – covers[i] indicates which indices j have alpha[j,:] participate in the i-th AGE cone. Automatically constructed in a presolve phase, if not provided.

  • c (Expression or None) – When provided, this DualSageCone instance will compile to a constraint to ensure that v is a valid dual variable to the constraint that \(c \in C_{\mathrm{SAGE}}(\alpha, X)\), If we have have information about the sign of a component of c, then it is possible to reduce the number of coniclifts primitives needed to represent this constraint.

alpha

The rows of alpha are the exponents defining this primal SAGE cone.

Type

ndarray

v

The vector subject to this dual SAGE cone constraint.

Type

Expression

X

If None, then this constraint represents a dual ordinary-SAGE cone.

Type

SigDomain

mu_vars

mu_vars[i] is the auxiliary variable associated with the i-th dual AGE cone. These variables are of shape mu_vars[i].size == alpha.shape[1]. The most basic solution recovery algorithm takes these variables, and considers points x of the form x = mu_vars[i].value / self.v[i].value.

Type

Dict[int, Variable]

ech

A simple wrapper around the constructor argument covers. Manages validation of covers when provided, and manages construction of covers when a user does not provide it.

Type

ExpCoverHelper

name

Uniquely identifies this Constraint in the model where it appears. Serves as a suffix for the name of any auxiliary Variable created when compiling to the coniclifts standard.

Type

str

The compiler interface

Up until now we have only described coniclifts as a tool for creating optimization problems. However, coniclifts’ more fundamental use is to exploit the following fact: for every convex set \(X \subset R^n\), there exists a matrix \(A \in R^{k \times m}\) , a vector \(b \in R^k\), and a convex cone \(K \subset R^k\) so that \(X = \{ (x_1,\ldots,x_n) \,:\, A x + b \in K, x \in R^m \}\). Coniclifts compiles all optimization problems into this standard form, where \(K\) is a product of elementary convex cones

  1. The zero cone.

  2. The nonnegative orthant.

  3. The exponential cone.

  4. The second-order cone.

  5. The vectorized positive semidefinite cone.

Crucially, coniclifts provides a means to map back and forth between models specified in high-level syntax, and models which exist in a flattened conic form using only primitives above.

The most important function in coniclifts’ compilation process is given below. The final return argument mentions “ScalarVariable” objects, which users of coniclifts need not interact with directly.

sageopt.coniclifts.compilers.compile_constrained_system(constraints)

Construct a flattened conic representation of the set of variable values which satisfy the constraints (“the feasible set”). Return the flattened representation of the feasible set, data structures for mapping the vectorized representation back to Variables, and a list of all Variables needed to represent the feasible set.

The final return argument (svid2col) can likely be ignored by end-users.

Parameters

constraints (list of coniclifts.Constraint) –

Returns

  • A (CSC-format sparse matrix) – The matrix appearing in the flattened representation of the feasible set.

  • b (ndarray) – The offset vector appearing in the flattened representation of the feasible set.

  • K (list of coniclifts Cone objects) – The cartesian product of these cones (in order) defines the convex cone appearing in the flattened representation of the feasible set.

  • variable_map (Dict[str, ndarray]) – A map from a Variable’s name field to a numpy array. If myvar is a coniclifts Variable appearing in the system defined by constraints, then a point x satisfying \(A x + b \in K\) maps to a feasible value for myvar by

    x0 = np.hstack([x, 0])
    myvar_val = x0[variable_map[myvar.name]]
    

    In particular, we guarantee myvar.shape == variable_map[myvar.name].shape. Augmenting x by zero to create x0 reflects a convention that if a component of a Variable does not affect the constraints, that component is automatically assigned the value zero.

  • variables (list of coniclifts.Variable) – All proper Variable objects appearing in the constraint set, including any auxiliary variables introduced to obtain a flattened conic system.

  • svid2col (Dict[int, int]) – A map from a ScalarVariable’s id to the index of the column in A where the ScalarVariable participates in the conic system. If the given ScalarVariable does not participate in the conic system, its id maps to -1.