6.4 Power Cone Optimization¶
The structure of a typical conic optimization problem is
(see Sec. 12 (Problem Formulation and Solutions) for detailed formulations). Here we discuss how to set-up problems with the primal/dual power cones.
MOSEK supports the primal and dual power cones, defined as below:
Primal power cone:
\[\POW_n^{\alpha_k} = \left\{ x\in \real^n ~:~ \prod_{i=0}^{n_\ell-1} x_i^{\beta_i} \geq \sqrt{\sum_{j=n_\ell}^{n-1}x_j^2},\ x_0\ldots,x_{n_\ell-1}\geq 0 \right\}\]where \(s = \sum_i \alpha_i\) and \(\beta_i = \alpha_i / s\), so that \(\sum_i \beta_i=1\).
Dual power cone:
\[(\POW_n^{\alpha_k}) = \left\{ x\in \real^n ~:~ \prod_{i=0}^{n_\ell-1} \left(\frac{x_i}{\beta_i}\right)^{\beta_i} \geq \sqrt{\sum_{j=n_\ell}^{n-1}x_j^2},\ x_0\ldots,x_{n_\ell-1}\geq 0 \right\}\]where \(s = \sum_i \alpha_i\) and \(\beta_i = \alpha_i / s\), so that \(\sum_i \beta_i=1\).
Perhaps the most important special case is the three-dimensional power cone family:
which has the corresponding dual cone:
For example, the conic constraint \((x,y,z)\in\POW_3^{0.25,0.75}\) is equivalent to \(x^{0.25}y^{0.75}\geq |z|\), or simply \(xy^3\geq z^4\) with \(x,y\geq 0\).
For other types of cones supported by MOSEK, see Sec. 15.10 (Supported domains) and the other tutorials in this chapter. Different cone types can appear together in one optimization problem.
6.4.1 Example POW1¶
Consider the following optimization problem which involves powers of variables:
We convert (6.12) into affine conic form using auxiliary variables as bounds for the power expressions:
The two conic constraints shown in (6.13) can be expressed in the ACC form as shown in (6.14):
(6.14)¶\[\begin{split}\left[\begin{array}{ccccc}1&0&0&0&0\\0&1&0&0&0\\0&0&0&1&0\\0&0&1&0&0\\0&0&0&0&0\\0&0&0&0&1\end{array}\right] \left[\begin{array}{c}x_0\\x_1\\x_2\\x_3\\x_4\end{array}\right] + \left[\begin{array}{c}0\\0\\0\\0\\1\\0\end{array}\right] \in \POW^{0.2,0.8}_3 \times \POW^{0.4,0.6}_3.\end{split}\]
Setting up the linear part
The linear parts (constraints, variables, objective) are set up using exactly the same methods as for linear problems, and we refer to Sec. 6.1 (Linear Optimization) for all the details. The same applies to technical aspects such as defining an optimization task, retrieving the solution and so on.
Setting up the conic constraints
In order to append the conic constraints we first input the matrix \(\afef\) and vector \(\afeg\) which together determine all the six affine expressions appearing in the conic constraints of (6.13)
/* Set the non-zero entries of the F matrix */
r = MSK_putafefentrylist(task, f_nnz, afeidx, varidx, f_val);
/* Set the non-zero element of the g vector */
if (r == MSK_RES_OK)
r = MSK_putafeg(task, 4, g);
/* Append the primal power cone domains with their respective parameter values. */
if (r == MSK_RES_OK)
r = MSK_appendprimalpowerconedomain(task, 3, 2, alpha_1, domidx);
if (r == MSK_RES_OK)
r = MSK_appendprimalpowerconedomain(task, 3, 2, alpha_2, domidx+1);
/* Append two ACCs made up of the AFEs and the domains defined above. */
if (r == MSK_RES_OK)
r = MSK_appendaccsseq(task, numacc, domidx, numafe, afeidx[0], NULL);
Following that, each of the affine conic constraints is appended using the function MSK_appendacc. The first argument selects the domain, which must be appended before being used, and must have the dimension matching the number of affine expressions appearing in the constraint. In the first case we append the power cone determined by the first three rows of \(\afef\) and \(\afeg\) while in the second call we use the remaining three rows of \(\afef\) and \(\afeg\).
Variants of this method are available to append multiple ACCs at a time. It is also possible to define the matrix \(\afef\) using a variety of methods (row after row, column by column, individual entries, etc.) similarly as for the linear constraint matrix \(A\).
For a more thorough exposition of the affine expression storage (AFE) matrix \(\afef\) and vector \(\afeg\) see Sec. 6.2 (From Linear to Conic Optimization).
Source code
#include <stdio.h>
#include "mosek.h" /* Include the MOSEK definition file. */
static void MSKAPI printstr(void *handle,
const char str[])
{
printf("%s", str);
} /* printstr */
int main(int argc, const char *argv[])
{
MSKrescodee r;
const MSKint32t numvar = 5,
numcon = 1;
const MSKint64t numafe = 6,
numacc = 2,
f_nnz = 5;
MSKboundkeye bkx[5];
double blx[5], bux[5];
double val[] = { 1.0, 1.0, -1.0 };
MSKint32t sub[] = { 3, 4, 0 };
double aval[] = { 1.0, 1.0, 0.5 };
MSKint32t asub[] = { 0, 1, 2 };
MSKint64t afeidx[] = {0, 1, 2, 3, 5};
MSKint32t varidx[] = {0, 1, 3, 2, 4};
MSKrealt f_val[] = {1, 1, 1, 1, 1},
g = 1;
const MSKrealt alpha_1[]= {0.2, 0.8},
alpha_2[]= {0.4, 0.6};
MSKint64t domidx[] = {0, 0};
MSKint32t i, j;
MSKenv_t env = NULL;
MSKtask_t task = NULL;
/* Create the mosek environment. */
r = MSK_makeenv(&env, NULL);
if (r == MSK_RES_OK)
{
/* Create the optimization task. */
r = MSK_maketask(env, numcon, numvar, &task);
if (r == MSK_RES_OK)
{
MSK_linkfunctotaskstream(task, MSK_STREAM_LOG, NULL, printstr);
/* Append 'numcon' empty constraints.
The constraints will initially have no bounds. */
if (r == MSK_RES_OK)
r = MSK_appendcons(task, numcon);
/* Append 'numvar' variables.
The variables will initially be fixed at zero (x=0). */
if (r == MSK_RES_OK)
r = MSK_appendvars(task, numvar);
/* Append 'numafe' affine expressions.
The affine expressions will initially be empty. */
if (r == MSK_RES_OK)
r = MSK_appendafes(task, numafe);
/* Set up the linear part */
MSK_putclist(task, 3, sub, val);
MSK_putarow(task, 0, 3, asub, aval);
MSK_putconbound(task, 0, MSK_BK_FX, 2.0, 2.0);
for(i=0;i<5;i++)
bkx[i] = MSK_BK_FR, blx[i] = -MSK_INFINITY, bux[i] = MSK_INFINITY;
MSK_putvarboundslice(task, 0, numvar, bkx, blx, bux);
if (r == MSK_RES_OK)
{
/* Set the non-zero entries of the F matrix */
r = MSK_putafefentrylist(task, f_nnz, afeidx, varidx, f_val);
/* Set the non-zero element of the g vector */
if (r == MSK_RES_OK)
r = MSK_putafeg(task, 4, g);
/* Append the primal power cone domains with their respective parameter values. */
if (r == MSK_RES_OK)
r = MSK_appendprimalpowerconedomain(task, 3, 2, alpha_1, domidx);
if (r == MSK_RES_OK)
r = MSK_appendprimalpowerconedomain(task, 3, 2, alpha_2, domidx+1);
/* Append two ACCs made up of the AFEs and the domains defined above. */
if (r == MSK_RES_OK)
r = MSK_appendaccsseq(task, numacc, domidx, numafe, afeidx[0], NULL);
}
MSK_putobjsense(task, MSK_OBJECTIVE_SENSE_MAXIMIZE);
if (r == MSK_RES_OK)
{
MSKrescodee trmcode;
/* Run optimizer */
r = MSK_optimizetrm(task, &trmcode);
/* Print a summary containing information
about the solution for debugging purposes*/
MSK_solutionsummary(task, MSK_STREAM_MSG);
if (r == MSK_RES_OK)
{
MSKsolstae solsta;
MSK_getsolsta(task, MSK_SOL_ITR, &solsta);
switch (solsta)
{
case MSK_SOL_STA_OPTIMAL:
{
double *xx = NULL;
xx = calloc(numvar, sizeof(double));
if (xx)
{
MSK_getxx(task,
MSK_SOL_ITR, /* Request the interior solution. */
xx);
printf("Optimal primal solution\n");
for (j = 0; j < 3; ++j)
printf("x[%d]: %e\n", j, xx[j]);
}
else
{
r = MSK_RES_ERR_SPACE;
}
free(xx);
}
break;
case MSK_SOL_STA_DUAL_INFEAS_CER:
case MSK_SOL_STA_PRIM_INFEAS_CER:
printf("Primal or dual infeasibility certificate found.\n");
break;
case MSK_SOL_STA_UNKNOWN:
printf("The status of the solution could not be determined. Termination code: %d.\n", trmcode);
break;
default:
printf("Other solution status.");
break;
}
}
else
{
printf("Error while optimizing.\n");
}
}
if (r != MSK_RES_OK)
{
/* In case of an error print error code and description. */
char symname[MSK_MAX_STR_LEN];
char desc[MSK_MAX_STR_LEN];
printf("An error occurred while optimizing.\n");
MSK_getcodedesc(r,
symname,
desc);
printf("Error %s - '%s'\n", symname, desc);
}
}
/* Delete the task and the associated data. */
MSK_deletetask(&task);
}
/* Delete the environment and the associated data. */
MSK_deleteenv(&env);
return (r);
} /* main */
