# Tutorial¶

First of all, let us import PICOS, CVXOPT, and a pretty printing helper

```
>>> import picos as pic
>>> import cvxopt as cvx
>>> from pprint import pprint
```

We now generate some arbitrary data, that we will use in this tutorial.

```
>>> pairs = [(0,2), (1,4), (1,3), (3,2), (0,4), (2,4)]
>>> b = ([0 ,2 ,0 ,3 ],
... [1 ,1 ,0 ,5 ],
... [-1,0 ,2 ,4 ],
... [0 ,0 ,-2,-1],
... [1 ,1 ,0 ,0 ])
>>> A = [] # A list of 2×4 matrices.
>>> for i in range(5):
... A.append(cvx.matrix(range(i-3,i+5),(2,4)))
>>> D = {'Peter': 12,
... 'Bob' : 4,
... 'Betty': 7,
... 'Elisa': 14}
```

Let us now create an instance `prob`

of an optimization problem

```
>>> prob = pic.Problem()
```

## Output settings¶

PICOS makes heavy use of unicode symbols to generate pretty output. If you find
that some of these symbols are not available on your terminal, you can call
`picos.ascii()`

or
`picos.latin1()`

to restrict the charset used.

```
>>> # Create a dummy symmetric matrix variable to showcase different charsets.
>>> X = pic.Problem().add_variable('X', (2,2), vtype="symmetric")
>>> # Show the default text representation of the SDP constraint X ≽ 0.
>>> print(X >> 0)
X ≽ 0
>>> # Restrict to ISO 8859-1 (aka Latin-1).
>>> pic.latin1()
>>> print(X >> 0)
X » 0
>>> # Restrict to pure ASCII.
>>> pic.ascii()
>>> print(X >> 0)
X >> 0
```

In addition to complete charset changes, PICOS allows you to fine-tune how certain expressions are formatted. See picos.glyphs for more.

```
>>> pic.glyphs.psdge.template = "{} - {} is positive semidefinite"
>>> print(X >> 0)
X - 0 is positive semidefinite
```

For the sake of this tutorial, we return to the full beauty of unicode.

```
>>> pic.default_charset() # Revert the changes; the same as pic.unicode().
```

## Variables¶

We will now create the variables of our optimization problem. This is done
by calling the method `add_variable()`

.
This function adds an instance of the class `Variable`

in the dictionary `prob.variables`

, and returns a reference
to the freshly added variable.
As we will next see, we
can use
this `Variable`

to form affine and quadratic expressions.

```
>>> t = prob.add_variable('t',1) #a scalar
>>> x = prob.add_variable('x',4) #a column vector
>>> Y = prob.add_variable('Y',(2,4)) #a matrix
>>> Z = []
>>> for i in range(5):
... Z.append(prob.add_variable('Z[{0}]'.format(i),(4,2))) # a list of 5 matrices
>>> w = {}
>>> for p in pairs: #a dictionary of (scalar) binary variables, indexed by our pairs
... w[p] = prob.add_variable('w[{0}]'.format(p),1 , vtype='binary')
```

Now, if we try to display a variable, here is what we get:

```
>>> w[2,4]
<1×1 Binary Variable: w[(2, 4)]>
>>> Y
<2×4 Continuous Variable: Y>
```

Also note the use of the
attributes `name`

, `value`

,
`size`

, and `vtype`

:

```
>>> w[2,4].vtype
'binary'
>>> x.vtype
'continuous'
>>> x.vtype='integer'
>>> x
<4×1 Integer Variable: x>
>>> x.size
(4, 1)
>>> Z[1].value = A[0].T
>>> Z[0].is_valued()
False
>>> Z[1].is_valued()
True
>>> Z[2].name
'Z[2]'
```

The admissible values for the `vtype`

attribute are documented in `add_variable()`

.

## Affine Expressions¶

We will now use our variables to create some affine expressions,
which are stored as instance of the class `AffinExp`

,
and will be the
core to define an optimization problem. Most python operators have been overloaded
to work with instances of `AffinExp`

(a list of available overloaded operators can be found in the
doc of `AffinExp`

). For example,
you can form the sum of two variables by writing:

```
>>> Z[0]+Z[3]
<4×2 Affine Expression: Z[0] + Z[3]>
```

The transposition of an affine expression is done by appending `.T`

:

```
>>> x
<4×1 Integer Variable: x>
>>> x.T
<1×4 Affine Expression: xᵀ>
```

### Parameters as constant affine expressions¶

It is also possible to form affine expressions by using parameters
stored in data structures such as a `list`

or a `cvxopt matrix`

(In fact, any type that is recognizable by the function `retrieve_matrix`

).

```
>>> x + b[0]
<4×1 Affine Expression: x + [4×1]>
>>> x.T + b[0]
<1×4 Affine Expression: xᵀ + [1×4]>
>>> A[0] * Z[0] + A[4] * Z[4]
<2×2 Affine Expression: [2×4]·Z[0] + [2×4]·Z[4]>
```

In the above example, you see that the list `b[0]`

was correctly converted into
a vector in the first expression, and into
a vector in the second one. This is because the overloaded
operators always try to convert the data into matrices of the appropriate size.

If you want to have better-looking string representations of your affine expressions,
you will need to convert the parameters into constant affine expressions. This can be done
thanks to the function `new_param()`

:

```
>>> A = pic.new_param('A',A) #this creates a list of constant affine expressions [A[0],...,A[4]]
>>> b = pic.new_param('b',b) #this creates a list of constant affine expressions [b[0],...,b[4]]
>>> D = pic.new_param('D',D) #this creates a dictionary of constant AffExpr, indexed by 'Peter', 'Bob', ...
>>> alpha = pic.new_param('alpha',12) #a scalar parameter
```

```
>>> alpha
<1×1 Affine Expression: alpha>
>>> D['Betty']
<1×1 Affine Expression: D[Betty]>
>>> pprint(b)
[<4×1 Affine Expression: b[0]>,
<4×1 Affine Expression: b[1]>,
<4×1 Affine Expression: b[2]>,
<4×1 Affine Expression: b[3]>,
<4×1 Affine Expression: b[4]>]
>>> print(b[0])
[ 0.00e+00]
[ 2.00e+00]
[ 0.00e+00]
[ 3.00e+00]
<BLANKLINE>
```

The above example also illustrates that when a *valued* affine expression
is printed, then its value is displayed instead of a symbolic description.
Note that the constant affine expressions, as `b[0]`

in the above example,
are always *valued*.
To assign a value to a non-constant `AffinExp`

,
you must set the `value`

property of
every variable involved in the affine expression.

```
>>> x_minus_1 = x - 1
>>> x_minus_1 #note that 1 was recognized as the (4x1)-vector with all ones
<4×1 Affine Expression: x - [1]>
>>> print(x_minus_1)
x - [1]
>>> x_minus_1.is_valued()
False
>>> x.value = [0,1,2,-1]
>>> x_minus_1.is_valued()
True
>>> print(x_minus_1)
[-1.00e+00]
[ 0.00e+00]
[ 1.00e+00]
[-2.00e+00]
<BLANKLINE>
```

We also point out that `new_param()`

converts lists into vectors and lists of lists into matrices (given
in row major order).
In contrast, tuples are converted into list of affine expressions:

```
>>> pic.new_param('vect',[1,2,3]) # [1,2,3] is converted into a vector of dimension 3
<3×1 Affine Expression: vect>
>>> pic.new_param('mat',[[1,2,3],[4,5,6]]) # [[1,2,3],[4,5,6]] is converted into a (2x3)-matrix
<2×3 Affine Expression: mat>
>>> pic.new_param('list_of_scalars',(1,2,3)) # (1,2,3) is converted into a list of 3 scalar parameters #doctest: +NORMALIZE_WHITESPACE
[<1×1 Affine Expression: list_of_scalars[0]>,
<1×1 Affine Expression: list_of_scalars[1]>,
<1×1 Affine Expression: list_of_scalars[2]>]
>>> pic.new_param('list_of_vectors',([1,2,3],[4,5,6])) # ([1,2,3],[4,5,6]) is converted into a list of 2 vector parameters #doctest: +NORMALIZE_WHITESPACE
[<3×1 Affine Expression: list_of_vectors[0]>,
<3×1 Affine Expression: list_of_vectors[1]>]
```

### Overloaded operators¶

OK, so now we have some variables (`t`

, `x`

, `w`

, `Y`

, and `Z`

)
and some parameters (`A`

, `b`

, `D`

and `alpha`

). Let us create some
affine expressions with them.

```
>>> A[0] * Z[0] #left multiplication
<2×2 Affine Expression: A[0]·Z[0]>
>>> Z[0] * A[0] #right multiplication
<4×4 Affine Expression: Z[0]·A[0]>
>>> A[1] * Z[0] * A[2] #left and right multiplication
<2×4 Affine Expression: A[1]·Z[0]·A[2]>
>>> alpha*Y #scalar multiplication
<2×4 Affine Expression: alpha·Y>
>>> t/b[1][3] - D['Bob'] #division by a scalar and substraction
<1×1 Affine Expression: t/b[1][3] - D[Bob]>
>>> ( b[2] | x ) #dot product
<1×1 Affine Expression: ⟨b[2], x⟩>
>>> ( A[3] | Y ) #generalized dot product for matrices: (A|B)=trace(A*B.T)
<1×1 Affine Expression: ⟨A[3], Y⟩>
>>> b[1]^x #hadamard (element-wise) product
<4×1 Affine Expression: b[1]⊙x>
```

We can also take some subelements of affine expressions, by using the standard syntax of python slices:

```
>>> b[1][1:3] #2d and 3rd elements of b[1]
<2×1 Affine Expression: b[1][1:3]>
>>> Y[1,:] #2d row of Y
<1×4 Affine Expression: Y[1,:]>
>>> x[-1] #last element of x
<1×1 Affine Expression: x[-1]>
>>> A[2][:,1:3]*Y[:,-2::-2] #extended slicing with (negative) steps is allowed
<2×2 Affine Expression: A[2][:,1:3]·Y[:,-2::-2]>
```

In the last example, we keep only the second and third columns of `A[2]`

, and
the columns of `Y`

with an even index, considered in the reverse order.
To concatenate affine expressions, the operators `//`

(vertical concatenation)
and `&`

(horizontal concatenation) have been overloaded:

```
>>> (b[1] & b[2] & x & A[0].T*A[0]*x) // x.T
<5×4 Affine Expression: [b[1], b[2], x, A[0]ᵀ·A[0]·x; xᵀ]>
```

When a scalar is added/substracted to a matrix or a vector, we interprete it as an elementwise addition of the scalar to every element of the matrix or vector.

```
>>> 5*x - alpha
<4×1 Affine Expression: 5·x - [alpha]>
```

### Summing Affine Expressions¶

You can take the advantage of python syntax to create sums of affine expressions:

```
>>> sum([A[i]*Z[i] for i in range(5)])
<2×2 Affine Expression: A[0]·Z[0] + A[1]·Z[1] + A[2]·Z[2] + A[3]·Z[3] + A[4]·Z[4]>
```

This works, but you might have very long string representations if there are a lot
of summands. So you’d better use the function `picos.sum()`

):

```
>>> pic.sum([A[i]*Z[i] for i in range(5)])
<2×2 Affine Expression: ∑(A[i]·Z[i] : i ∈ [0…4])>
```

## Objective function¶

The objective function of the problem
can be defined with the function
`set_objective()`

.
Its first argument should be `'max'`

, `'min'`

or
`'find'`

(*for feasibility problems*),
and the second argument should be a scalar expression:

```
>>> prob.set_objective('max',( A[0] | Y )-t)
>>> print(prob) #doctest: +NORMALIZE_WHITESPACE
---------------------
optimization problem (Mixed-Integer Unconstrained Problem):
59 variables, 0 affine constraints
<BLANKLINE>
Y : (2, 4), continuous
Z : list of 5 variables, (4, 2), continuous
t : (1, 1), continuous
w : dict of 6 variables, (1, 1), binary
x : (4, 1), integer
<BLANKLINE>
maximize ⟨A[0], Y⟩ - t
such that
[]
---------------------
```

With this example, you see what happens when a problem is printed: the list of optimization variables is displayed, then the objective function and finally a list of constraints (in the case above, there is no constraint).

## Norm of an affine Expression¶

The norm of an affine expression is an overload of the `abs()`

function.
If `x`

is an affine expression, `abs(x)`

is its Euclidean norm .

```
>>> abs(x)
<Norm of a 4×1 Expression: ‖x‖>
```

In the case where the affine expression is a matrix, `abs()`

returns its
Frobenius norm, defined as .

```
>>> abs(Z[1]-2*A[0].T)
<Norm of a 4×2 Expression: ‖Z[1] - 2·A[0]ᵀ‖>
```

Note that the absolute value of a scalar expression is stored as a norm:

```
>>> abs(t)
<Norm of a 1×1 Expression: ‖t‖>
```

However, a scalar constraint of the form is handled as two linear constraints by PICOS, and so a problem with the latter constraint can be solved even if you do not have a SOCP solver available.

It is also possible to use other norms in picos, cf. this paragraph .

## Quadratic Expressions¶

Quadratic expressions can be formed in several ways:

```
>>> t**2 - x[1]*x[2] + 2*t - alpha #sum of linear and quadratic terms
<Quadratic Expression: t² - x[1]·x[2] + 2·t - alpha>
>>> (x[1]-2) * (t+4) #product of two affine expressions
<Quadratic Expression: (x[1] - 2)·(t + 4)>
>>> Y[0,:]*x #Row vector multiplied by column vector
<Quadratic Expression: Y[0,:]·x>
>>> (x +2 | Z[1][:,1]) #scalar product of affine expressions
<Quadratic Expression: ⟨x + [2], Z[1][:,1]⟩>
>>> abs(x)**2 #recall that abs(x) is the euclidean norm of x
<Quadratic Expression: ‖x‖²>
>>> (t & alpha) * A[1] * x #quadratic form
<Quadratic Expression: [t, alpha]·A[1]·x>
```

It is not possible (yet) to make a multidimensional quadratic expression.

## Constraints¶

A constraint takes the form of two expressions separated by a relation operator.

### Linear (in)equalities¶

Linear (in)equalities are understood elementwise. **The strict operators**
`<`

**and** `>`

**denote weak inequalities** (*less or equal than*
and *larger or equal than*). For example:

```
>>> (1|x) < 2 #sum of the x[i] less or equal than 2
<1×1 Affine Constraint: ⟨[1], x⟩ ≤ 2>
>>> Z[0] * A[0] > b[1]*b[2].T #A 4x4-elementwise inequality
<4×4 Affine Constraint: Z[0]·A[0] ≥ b[1]·b[2]ᵀ>
>>> pic.sum([A[i]*Z[i] for i in range(5)]) == 0 #A 2x2 equality. The RHS is the all-zero matrix
<2×2 Affine Constraint: ∑(A[i]·Z[i] : i ∈ [0…4]) = 0>
```

Constraints can be added in the problem with the function
`add_constraint()`

:

```
>>> constraints = []
>>> for i in range(1,5):
... constraints.append(prob.add_constraint(Z[i]==Z[i-1]+Y.T))
>>> print(prob) #doctest: +NORMALIZE_WHITESPACE
---------------------
optimization problem (MILP):
59 variables, 32 affine constraints
<BLANKLINE>
Y : (2, 4), continuous
Z : list of 5 variables, (4, 2), continuous
t : (1, 1), continuous
w : dict of 6 variables, (1, 1), binary
x : (4, 1), integer
<BLANKLINE>
maximize ⟨A[0], Y⟩ - t
such that
Z[1] = Z[0] + Yᵀ
Z[2] = Z[1] + Yᵀ
Z[3] = Z[2] + Yᵀ
Z[4] = Z[3] + Yᵀ
---------------------
```

The constraints of the problem can be accessed via the reference returned by
`add_constraint()`

or with the function
`get_constraint()`

:

```
>>> prob.get_constraint(2)
<4×2 Affine Constraint: Z[3] = Z[2] + Yᵀ>
>>> constraints[2] is prob.get_constraint(2)
True
```

### Groupping constraints¶

In order to have a more compact string representation of the problem, it is
advised to use the function `add_list_of_constraints()`

, which works similarly as the
function `sum()`

.

```
>>> # We first remove all constraints previously added.
>>> prob.remove_all_constraints()
>>> # Then we add a single constraint.
>>> C = prob.add_constraint(Y > 0)
>>> # And then the same list of constraints as above.
>>> CL = prob.add_list_of_constraints(
... [Z[i] == Z[i-1] + Y.T for i in range(1,5)])
>>> print(prob) #doctest: +NORMALIZE_WHITESPACE
---------------------
optimization problem (MILP):
59 variables, 40 affine constraints
<BLANKLINE>
Y : (2, 4), continuous
Z : list of 5 variables, (4, 2), continuous
t : (1, 1), continuous
w : dict of 6 variables, (1, 1), binary
x : (4, 1), integer
<BLANKLINE>
maximize ⟨A[0], Y⟩ - t
such that
Y ≥ 0
Z[i+1] = Z[i] + Yᵀ ∀ i ∈ [0…3]
---------------------
```

Now, the constraint `Z[3] == Z[2] + Y.T`

, which has been entered
in 4th position, can either be accessed by `prob.get_constraint(3)`

(`3`

because constraints are numbered from `0`

), or by

```
>>> prob.get_constraint((1,2))
<4×2 Affine Constraint: Z[3] = Z[2] + Yᵀ>
```

where `(1,2)`

means *the 3rd constraint of the 2d group of constraints*,
with zero-based numbering.

Similarly, the constraint `Y ≥ 0`

can be accessed by
`prob.get_constraint(0)`

(first constraint),
`prob.get_constraint((0,0))`

(first constraint of the first group), or
`prob.get_constraint((0,))`

(unique constraint of the first group).

### Removing constraints¶

It can be useful to remove some constraints, especially for
dynamic approaches such as column generation. Re-creating an instance from
scratch after each iteration would be inefficient. Instead, PICOS allows
one to modify the solver instance and to re-solve it on the fly
(for `mosek`

, `cplex`

and `gurobi`

). It is not possible to change
directly a constraint, but you can delete a constraint from model, and
then re-add the modified constraint.

The next code shows an example in which a variable `x2`

is added to the model,
which appears as `+3*x2`

in the objective function and as `+x2`

in the LHS of a constraint.

```
>>> prb = pic.Problem()
>>> x1 = prb.add_variable('x1',1)
>>> lhs = 2*x1
>>> obj = 5*x1
>>> cons = prb.add_constraint(lhs <= 1)
>>> sol = prb.maximize(obj,verbose=0,solver='cvxopt')
>>> #--------------------------------------
>>> #at this place, the user can use his favorite method to solve the 'pricing problem'.
>>> #Let us assume that this phase suggests to add a new variable 'x2' in the model,
>>> #which appears as `+3*x2` in the objective function, and as `+x2` in the LHS of the constraint.
>>> x2 = prb.add_variable('x2',1)
>>> lhs += x2
>>> obj += 3*x2
>>> cons.delete()
>>> newcons = prb.add_constraint(lhs <= 1)
>>> print(prb) #doctest: +NORMALIZE_WHITESPACE
---------------------
optimization problem (LP):
2 variables, 1 affine constraints
<BLANKLINE>
x1 : (1, 1), continuous
x2 : (1, 1), continuous
<BLANKLINE>
maximize 5·x1 + 3·x2
such that
2·x1 + x2 ≤ 1
---------------------
```

### Flow constraints in Graphs¶

Flow constraints in graphs are entered using a Networkx Graph. The following example finds a (trivial) maximal flow from `'S'`

to `'T'`

in `G`

.

```
>>> import networkx as nx
>>> G = nx.DiGraph()
>>> G.add_edge('S','A', capacity=1)
>>> G.add_edge('A','B', capacity=1)
>>> G.add_edge('B','T', capacity=1)
>>> pb = pic.Problem()
>>> # Adding the flow variables
>>> f={}
>>> for e in G.edges():
... f[e]=pb.add_variable('f[{0}]'.format(e),1)
>>> # A variable for the value of the flow
>>> F = pb.add_variable('F',1)
>>> # Creating the flow constraint
>>> flowCons = pb.add_constraint(pic.flow_Constraint(
... G, f, 'S', 'T', F, capacity='capacity', graphName='G'))
>>> pb.set_objective('max',F)
>>> sol = pb.solve(verbose=0,solver='cvxopt')
>>> flow = pic.tools.eval_dict(f)
```

Picos allows you to define single source - multiple sinks problems. You can use the same syntax as for a single source - single sink problems. Just add a list of sinks and a list of flows instead.

```
import picos as pic
import networkx as nx
G=nx.DiGraph()
G.add_edge('S','A', capacity=2); G.add_edge('S','B', capacity=2)
G.add_edge('A','T1', capacity=2); G.add_edge('B','T2', capacity=2)
pbMultipleSinks=pic.Problem()
# Flow variable
f={}
for e in G.edges():
f[e]=pbMultipleSinks.add_variable('f[{0}]'.format(e),1)
# Flow value
F1=pbMultipleSinks.add_variable('F1',1)
F2=pbMultipleSinks.add_variable('F2',1)
flowCons = pic.flow_Constraint(
G, f, source='S', sink=['T1','T2'], capacity='capacity',
flow_value=[F1, F2], graphName='G')
pbMultipleSinks.add_constraint(flowCons)
pbMultipleSinks.set_objective('max',F1+F2)
# Solve the problem
pbMultipleSinks.solve(verbose=0,solver='cvxopt')
print(pbMultipleSinks)
print('The optimal flow F1 has value {:.1f}'.format(F1.value))
print('The optimal flow F2 has value {:.1f}'.format(F2.value))
```

```
---------------------
optimization problem (LP):
6 variables, 12 affine constraints
F1 : (1, 1), continuous
F2 : (1, 1), continuous
f : dict of 4 variables, (1, 1), continuous
maximize F1 + F2
such that
Feasible S-(T1, T2)-flow in G has value (F1, F2).
---------------------
The optimal flow F1 has value 2.0
The optimal flow F2 has value 2.0
```

A similar syntax can be used for multiple sources-single sink flows.

### Quadratic constraints¶

Quadratic inequalities are entered in the following way:

```
>>> t**2 > 2*t - alpha + x[1]*x[2]
<Quadratic Constraint: -t² + 2·t - alpha + x[1]·x[2] ≤ 0>
>>> (t & alpha) * A[1] * x + (x +2 | Z[1][:,1]) < 3*(1|Y)-alpha
<Quadratic Constraint: [t, alpha]·A[1]·x + ⟨x + [2], Z[1][:,1]⟩ - (3·⟨[1], Y⟩ - alpha) ≤ 0>
```

Note that PICOS does not check the convexity of convex constraints. The solver will raise an Exception if it does not support non-convex quadratics.

### Second Order Cone Constraints¶

There are two types of second order cone constraints supported in PICOS.

The constraints of the type , where is a scalar affine expression and is a multidimensional affine expression (possibly a matrix, in which case the norm is Frobenius). This inequality forces the vector to belong to a Lorrentz-Cone (also called

ice-cream cone).The constraints of the type , where and are scalar affine expressions and is a multidimensional affine expression, which constrain the vector inside a rotated version of the Lorretz cone.

A few examples:

```
>>> # A simple ice-cream cone constraint
>>> abs(x) < (2|x-1)
<5×1 SOC Constraint: ‖x‖ ≤ ⟨[2], x - [1]⟩>
>>> # SOC constraint with Frobenius norm
>>> abs(Y+Z[0].T) < t+alpha
<9×1 SOC Constraint: ‖Y + Z[0]ᵀ‖ ≤ t + alpha>
>>> # Rotated SOC constraint
>>> abs(Z[1][:,0])**2 < (2*t-alpha)*(x[2]-x[-1])
<6×1 RSOC Constraint: ‖Z[1][:,0]‖² ≤ (2·t - alpha)·(x[2] - x[-1])>
>>> # t**2 is understood as the squared norm of [t]
>>> t**2 < D['Elisa']+t
<3×1 RSOC Constraint: ‖t‖² ≤ D[Elisa] + t>
>>> # 1 is understood as the squared norm of [1]
>>> 1 < (t-1)*(x[2]+x[3])
<3×1 RSOC Constraint: ‖1^(1/2)‖² ≤ (t - 1)·(x[2] + x[3])>
```

### Semidefinite Constraints¶

Linear matrix inequalities (LMI) can be entered thanks to an overload of the operators
`<<`

and `>>`

. For example, the LMI

where is used to denote the Löwner ordering, is passed to PICOS by writing:

```
>>> pic.sum([x[i]*b[i]*b[i].T for i in range(4)]) >> b[4]*b[4].T
<4×4 LMI Constraint: ∑(x[i]·b[i]·b[i]ᵀ : i ∈ [0…3]) ≽ b[4]·b[4]ᵀ>
```

Note the difference with

```
>>> pic.sum([x[i]*b[i]*b[i].T for i in range(4)]) > b[4]*b[4].T
<4×4 Affine Constraint: ∑(x[i]·b[i]·b[i]ᵀ : i ∈ [0…3]) ≥ b[4]·b[4]ᵀ>
```

which yields an elementwise inequality.

For convenience, it is possible to add a symmetric matrix variable `X`

,
by specifying the option `vtype=symmetric`

. This has the effect to
store all the affine expressions which depend on `X`

as a function
of its lower triangular elements only.

```
>>> sdp = pic.Problem()
>>> X = sdp.add_variable('X',(4,4),vtype='symmetric')
>>> C = sdp.add_constraint(X >> 0)
>>> print(sdp) #doctest: +NORMALIZE_WHITESPACE
---------------------
optimization problem (SDP):
10 variables, 0 affine constraints, 10 vars in 1 SD cones
<BLANKLINE>
X : (4, 4), symmetric
<BLANKLINE>
find vars
such that
X ≽ 0
---------------------
```

In this example, you see indeed that the problem has 10=(4*5)/2 variables,
which correspond to the lower triangular elements of `X`

.

Warning

When a constraint of the form `A >> B`

is passed to PICOS, it is not
assumed that A-B is symmetric. Instead, the symmetric matrix whose lower
triangular elements are those of A-B is forced to be positive semidefnite.
So, in the cases where A-B is not implicitely forced to be symmetric, you
should add a constraint of the form `A-B==(A-B).T`

in the problem.

### Inequalities involving geometric means¶

It is possible to enter an inequality of the form

in PICOS, where is a scalar affine expression and is an affine expression of dimension (possibly a matrix, in which case is counted in column major order). This inequality is internally converted to an equivalent set of second order cone inequalities, by using standard techniques (cf. e.g. [1]).

Many convex constraints can be formulated using inequalities that involve
a geometric mean. For example, is equivalent
to , which can be entered in PICOS
thanks to the function `picos.geomean()`

:

```
>>> t < pic.geomean(t //x[1] //x[1] //1)
<Geometric Mean Constraint: t ≤ geomean([t; x[1]; x[1]; 1])>
```

Note that the latter example can also be passed to picos in a more simple way,
thanks to an overloading of the `**`

exponentiation operator:

```
>>> t < x[1]**(2./3)
<Power Constraint: x[1]^(2/3) ≥ t>
```

Inequalities involving geometric means are stored in a temporary object
of the class `GeoMeanConstraint`

,
which can be passed to a problem with
`add_constraint()`

:

```
>>> geom_ineq = prob.add_constraint(t<pic.geomean(Y[:6]))
```

The constraint has an attribute `constraints`

which contains auxiliary SOC
inequalities that are used internally to represent the geometric mean:

```
>>> pprint(geom_ineq.constraints)
(<3×1 RSOC Constraint: ‖u[1:0-1]‖² ≤ Y[:6][0]·Y[:6][1]>,
<3×1 RSOC Constraint: ‖u[1:2-3]‖² ≤ Y[:6][2]·Y[:6][3]>,
<3×1 RSOC Constraint: ‖u[1:4-5]‖² ≤ Y[:6][4]·Y[:6][5]>,
<3×1 RSOC Constraint: ‖u[2:0-3]‖² ≤ u[1:0-1]·u[1:2-3]>,
<3×1 RSOC Constraint: ‖u[2:4-x]‖² ≤ u[1:4-5]·t>,
<3×1 RSOC Constraint: ‖t‖² ≤ u[2:0-3]·u[2:4-x]>)
```

Likewise, auxiliary variables can be accessed via the `variables`

attribute:

```
>>> pprint(geom_ineq.variables)
(<1×1 Continuous Variable: u[1:0-1]>,
<1×1 Continuous Variable: u[1:2-3]>,
<1×1 Continuous Variable: u[1:4-5]>,
<1×1 Continuous Variable: u[2:0-3]>,
<1×1 Continuous Variable: u[2:4-x]>)
```

### Inequalities involving real powers or trace of matrix powers¶

As mentionned above, the `**`

exponentiation operator has been overloaded
to support real exponents. A rational approximation of the exponent is used,
and the inequality are internally reformulated as a set of equivalent SOC inequalities.
Note that only inequalities defining a convex regions can be passed:

```
>>> t**0.6666 > x[0]
<Power Constraint: t^(2/3) ≥ x[0]>
>>> t**-0.5 < x[0]
<Power Constraint: t^(-1/2) ≤ x[0]>
>>> try:
... t**-0.5 > x[0]
... except Exception as ex:
... print('Exception: '+str(ex)) #doctest: +NORMALIZE_WHITESPACE
Exception: May only lower bound a concave power (p >= 0 and p <= 1).
>>> t**2 < x[1]+x[2]
<3×1 RSOC Constraint: ‖t‖² ≤ x[1] + x[2]>
```

More generally, inequalities involving trace of matrix powers can be passed to PICOS,
by using the `picos.tracepow()`

function. The following example
creates the constraint

```
>>> pic.tracepow(x[0] * A[0]*A[0].T + x[2] * A[2]*A[2].T, 2.5) <= 3
<Trace of Power Constraint: trace((x[0]·A[0]·A[0]ᵀ + x[2]·A[2]·A[2]ᵀ)^(5/2)) ≤ 3>
```

Warning

when a power expression (resp. the trace of matrix power ) is used, the base is forced to be nonnegative (resp. the base is forced to be positive semidefinite) by picos.

When the exponent is , it is also possible to represent constraints of the form with SDPs, where , see [2].

```
>>> pic.tracepow(X, 0.6666, coef = A[0].T*A[0]) >= t
<Trace of Power Constraint: trace(A[0]ᵀ·A[0]·X^(2/3)) ≥ t>
```

As for geometric means, inequalities involving real powers yield their internal
representation via the `constraints`

and `variables`

attributes.

### Inequalities involving generalized p-norm¶

Inequalities of the form can be entered by using the
function `picos.norm()`

. This function is also defined for
by the usual formula .
The norm function is convex over for all , and
concave over the set of vectors with nonnegative coordinates for .

```
>>> pic.norm(x,3) < t
<p-Norm Constraint: ‖x‖_3 ≤ t>
>>> pic.norm(x,'inf') < 2
<p-Norm Constraint: ‖x‖_inf ≤ 2>
>>> pic.norm(x,0.5) > x[0]-x[1]
<Generalized p-Norm Constraint: ‖x‖_(1/2) ≥ x[0] - x[1]>
```

Warning

Note that when a constraint of the form `norm(x,p) >= t`

is entered (with ),
PICOS forces the vector `x`

to be nonnegative (componentwise).

Inequalities involving the generalized norm of
a matrix can also be handled with picos, cf. the documentation of
`picos.norm()`

.

As for geometric means, inequalities involving p-norms yield their internal
representation via the `constraints`

and `variables`

attributes.

### Inequalities involving the nth root of a determinant¶

The function `picos.detrootn()`

can be used to enter the -th root of the determinant of a
-symmetric positive semidefinite matrix:

```
>>> M = sdp.add_variable('M',(5,5),'symmetric')
>>> t < pic.detrootn(M)
<n-th Root of a Determinant Constraint: det(M)^(1/5) ≥ t>
```

Warning

Note that when a constraint of the form `t < pic.detrootn(M)`

is entered (with ),
PICOS forces the matrix `M`

to be positive semidefinite.

As for geometric means, inequalities involving the nth root of a determinant
yield their internal representation via the `constraints`

and `variables`

attributes.

## Set membership¶

Since Picos 1.0.2, there is a `Set`

class
that can be used to pass constraints as membership of an affine expression to a set.

Following sets are currently supported:

balls representing the set can be constructed with the function

`pic.ball()`

The standard simplex (scaled by a factor ) can be constructed with the function

`pic.simplex()`

Truncated simplexes and symmetrized Truncated simplexes can be constructed with the function

`pic.truncated_simplex()`

Membership of an affine expression to a set can be expressed with the overloaded operator `<<`

.
This returns a temporary object that can be passed to a picos problem with the function
`add_constraint()`

.

```
>>> x << pic.simplex(1)
<Standard Simplex Constraint: x ∈ {x ≥ 0 : ∑(x) ≤ 1}>
>>> x << pic.truncated_simplex(2)
<Truncated Simplex Constraint: x ∈ {0 ≤ x ≤ 1 : ∑(x) ≤ 2}>
>>> x << pic.truncated_simplex(2,sym=True)
<Symmetrized Truncated Simplex Constraint: x ∈ {-1 ≤ x ≤ 1 : ∑(|x|) ≤ 2}>
>>> x << pic.ball(3)
<5×1 SOC Constraint: ‖x‖ ≤ 3>
>>> pic.ball(2,'inf') >> x
<p-Norm Constraint: ‖x‖_inf ≤ 2>
>>> x << pic.ball(4,1.5)
<p-Norm Constraint: ‖x‖_(3/2) ≤ 4>
```

## Write a Problem to a file¶

It is possible to write a problem to a file, thanks to the
function `write_to_file()`

.
Several file formats and file writers are available, have a look at the doc
of `write_to_file()`

for more explanations.

Below is a *hello world* example, which writes a simple LP to a **.lp** file:

```
import picos as pic
prob = pic.Problem(verbose = False)
y = prob.add_variable('y',1)
x = prob.add_variable('x',1)
prob.add_constraint(x>1.5)
prob.add_constraint(y-x>0.7)
prob.set_objective('min',y)
#let first picos display the problem
print(prob)
print()
#now write the problem to a .lp file...
prob.write_to_file('helloworld.lp')
print()
#and display the content of the freshly created file:
print(open('helloworld.lp').read())
```

Generated output:

```
---------------------
optimization problem (LP):
2 variables, 2 affine constraints
x : (1, 1), continuous
y : (1, 1), continuous
minimize y
such that
x ≥ 1.5
y - x ≥ 0.7
---------------------
writing problem in helloworld.lp...
done.
\* file helloworld.lp generated by picos*\
Minimize
obj : 1 y
Subject To
in0 : -1 x <= -1.5
in1 : -1 y+ 1 x <= -0.7
Bounds
y free
x free
Generals
Binaries
End
```

## Solve a Problem¶

To solve a problem, you have to use the method `solve()`

of the class `Problem`

.
Alternatively, the functions `maximize(obj)`

and `minimize(obj)`

can be used to specify the objective function and call the solver in a single statement.
These method accept several options. In particular the solver can be specified by passing
an option of the form `solver='solver_name'`

. For a list of available
parameters with their default values, see the doc of the function
`set_all_options_to_default()`

.

Once a problem has been solved, the optimal values of the variables are
accessible with the `value`

property.
Depending on the solver, you
can also obtain the slack and the optimal dual variables
of the constraints thanks to the properties
`dual`

and
`slack`

of the class
`Constraint`

.
See the doc of `dual`

for more explanations
on the dual variables for second order cone programs (SOCP) and
semidefinite programs (SDP).

The class `Problem`

also has
two interesting properties: `type`

, which
indicates the class of the optimization problem (‘LP’, ‘SOCP’, ‘MIP’, ‘SDP’,…),
and `status`

, which indicates if the
problem has been solved (the default is `'unsolved'`

; after a call to
`solve()`

this property can take the value of any
code returned by a solver, such as `'optimal'`

, `'unbounded'`

, `'near-optimal'`

,
`'primal infeasible'`

, `'unknown'`

, …).

Below is a simple example, to solve the linear programm:

More examples can be found here.

```
P = pic.Problem()
A = pic.new_param('A', cvx.matrix([[1,1],[0,1]]) )
x = P.add_variable('x',2)
P.add_constraint(x[0]>x[1])
P.add_constraint(A*x<[3,4])
objective = 0.5 * x[0] + x[1]
P.set_objective('max', objective) #or directly P.maximize(objective)
#display the problem and solve it
print(P)
print('type: '+P.type)
print('status: '+P.status)
P.solve(verbose=0,solver='cvxopt')
print('status: '+P.status)
print()
#--------------------#
# objective value #
#--------------------#
print('The optimal value of this problem is:')
print(P.obj_value()) # "print(objective)" would also work, because objective is valued
print()
#--------------------#
# optimal variable #
#--------------------#
x_opt = x.value
print('The solution of the problem is:')
print(x_opt) # "print(x)" would also work, since x is now valued
#--------------------#
# slacks and duals #
#--------------------#
c0=P.get_constraint(0)
print('The dual of the constraint {} is:'.format(c0))
print(c0.dual)
print()
print('And its slack is:')
print(c0.slack)
print()
c1=P.get_constraint(1)
print('The dual of the constraint {} is:'.format(c1))
print(c1.dual)
print('And its slack is:')
print(c1.slack)
```

```
---------------------
optimization problem (LP):
2 variables, 3 affine constraints
x : (2, 1), continuous
maximize 0.5·x[0] + x[1]
such that
x[0] ≥ x[1]
A·x ≤ [2×1]
---------------------
type: LP
status: unsolved
status: optimal
The optimal value of this problem is:
3.0...
The solution of the problem is:
[ 2.00e+00]
[ 2.00e+00]
The dual of the constraint x[0] ≥ x[1] is:
0.25...
And its slack is:
1.8...
The dual of the constraint A·x ≤ [2×1] is:
[ 4.56e-10]
[ 7.50e-01]
And its slack is:
[ 1.00e+00]
[-8.71e-10]
```

## A note on dual variables¶

For second order cone constraints of the form , where is a vector of dimension , the dual variable is a vector of dimension of the form , where the -dimensional vector satisfies .

For rotated second order cone constraints of the form where is a vector of dimension , the dual variable is a vector of dimension of the form with , where the -dimensional vector satisfies .

In general, a linear problem with second order cone constraints (both standard and rotated) and semidefinite constraints can be written under the form:

where

are vectors of dimension ;

are scalars;

are vectors of dimension and are matrices of size ;

are vectors of dimension and are matrices of size ;

is a vector of dimension and is a matrix of size ;

is a vector of dimension and is a matrix of size ;

are symmetric matrices ().

Its dual problem can be written as:

where stands for the vector of dimension with on the th coordinate, and the dual variables are

When quering the dual of a constraint of the above primal problem, **picos will
return**

for the constraint ;

for the constraint ;

The -dimensional vector for the constraint

;

The -dimensional vector for the constraint

;

The symmetric positive definite matrix for the constraint

.

## References¶

“Applications of second-order cone programming”, M.S. Lobo, L. Vandenberghe, S. Boyd and H. Lebret,

Linear Algebra and its Applications, 284, p.193-228, 1998.“On the semidefinite representations of real functions applied to symmetric matrices”, G. Sagnol,

Linear Algebra and its Applications, 439(10), p.2829-2843, 2013.