/* Copyright 2024, Gurobi Optimization, LLC
This example formulates and solves the following simple model
with PWL constraints:
maximize
sum c[j] * x[j]
subject to
sum A[i,j] * x[j] <= 0, for i = 0, ..., m-1
sum y[j] <= 3
y[j] = pwl(x[j]), for j = 0, ..., n-1
x[j] free, y[j] >= 0, for j = 0, ..., n-1
where pwl(x) = 0, if x = 0
= 1+|x|, if x != 0
Note
1. sum pwl(x[j]) <= b is to bound x vector and also to favor sparse x vector.
Here b = 3 means that at most two x[j] can be nonzero and if two, then
sum x[j] <= 1
2. pwl(x) jumps from 1 to 0 and from 0 to 1, if x moves from negative 0 to 0,
then to positive 0, so we need three points at x = 0. x has infinite bounds
on both sides, the piece defined with two points (-1, 2) and (0, 1) can
extend x to -infinite. Overall we can use five points (-1, 2), (0, 1),
(0, 0), (0, 1) and (1, 2) to define y = pwl(x)
*/
#include <stdlib.h>
#include <stdio.h>
#include "gurobi_c.h"
int
main(int argc,
char *argv[])
{
GRBenv *env = NULL;
GRBmodel *model = NULL;
int *cbeg = NULL;
int *clen = NULL;
int *cind = NULL;
double *cval = NULL;
double *rhs = NULL;
char *sense = NULL;
double *lb = NULL;
double *obj = NULL;
int nz, i, j;
int status;
double objval;
int error = 0;
int n = 5;
int m = 5;
double c[] = { 0.5, 0.8, 0.5, 0.1, -1 };
double A[][5] = { {0, 0, 0, 1, -1},
{0, 0, 1, 1, -1},
{1, 1, 0, 0, -1},
{1, 0, 1, 0, -1},
{1, 0, 0, 1, -1} };
int npts = 5;
double xpts[] = {-1, 0, 0, 0, 1};
double ypts[] = {2, 1, 0, 1, 2};
/* Create environment */
error = GRBloadenv(&env, NULL);
if (error) goto QUIT;
/* Allocate memory and build the model */
nz = n; /* count nonzeros for n y variables */
for (i = 0; i < m; i++) {
for (j = 0; j < n; j++) {
if (A[i][j] != 0.0) nz++;
}
}
cbeg = (int *) malloc(2*n*sizeof(int));
clen = (int *) malloc(2*n*sizeof(int));
cind = (int *) malloc(nz*sizeof(int));
cval = (double *) malloc(nz*sizeof(double));
rhs = (double *) malloc((m+1)*sizeof(double));
sense = (char *) malloc((m+1)*sizeof(char));
lb = (double *) malloc(2*n*sizeof(double));
obj = (double *) malloc(2*n*sizeof(double));
for (j = 0; j < n; j++) {
/* for x variables */
lb[j] = -GRB_INFINITY;
obj[j] = c[j];
/* for y variables */
lb[j+n] = 0.0;
obj[j+n] = 0.0;
}
for (i = 0; i < m; i++) {
rhs[i] = 0.0;
sense[i] = GRB_LESS_EQUAL;
}
sense[m] = GRB_LESS_EQUAL;
rhs[m] = 3;
nz = 0;
for (j = 0; j < n; j++) {
cbeg[j] = nz;
for (i = 0; i < m; i++) {
if (A[i][j] != 0.0 ) {
cind[nz] = i;
cval[nz] = A[i][j];
nz++;
}
}
clen[j] = nz - cbeg[j];
}
for (j = 0; j < n; j++) {
cbeg[n+j] = nz;
clen[n+j] = 1;
cind[nz] = m;
cval[nz] = 1.0;
nz++;
}
error = GRBloadmodel(env, &model, "gc_pwl_c", 2*n, m+1,
GRB_MAXIMIZE, 0.0, obj, sense, rhs,
cbeg, clen, cind, cval, lb, NULL,
NULL, NULL, NULL);
if (error) goto QUIT;
/* Add piecewise constraints */
for (j = 0; j < n; j++) {
error = GRBaddgenconstrPWL(model, NULL, j, n+j, npts, xpts, ypts);
if (error) goto QUIT;
}
/* Optimize model */
error = GRBoptimize(model);
if (error) goto QUIT;
for (j = 0; j < n; j++) {
double x;
error = GRBgetdblattrelement(model, "X", j, &x);
if (error) goto QUIT;
printf("x[%d] = %g\n", j, x);
}
/* Report the result */
error = GRBgetintattr(model, GRB_INT_ATTR_STATUS, &status);
if (error) goto QUIT;
if (status != GRB_OPTIMAL) {
fprintf(stderr, "Error: it isn't optimal\n");
goto QUIT;
}
error = GRBgetdblattr(model, GRB_DBL_ATTR_OBJVAL, &objval);
if (error) goto QUIT;
printf("Obj: %g\n", objval);
QUIT:
/* Error reporting */
if (error) {
printf("ERROR: %s\n", GRBgeterrormsg(env));
exit(1);
}
/* Free data */
if (cbeg) free(cbeg);
if (clen) free(clen);
if (cind) free(cind);
if (cval) free(cval);
if (rhs) free(rhs);
if (sense) free(sense);
if (lb) free(lb);
if (obj) free(obj);
/* Free model */
GRBfreemodel(model);
/* Free environment */
GRBfreeenv(env);
return 0;
}