/***********************************************************************
Xpress Optimizer Examples
=========================
file Polygon_userfunc.c
```````````````````````
Implement the polygon example using a user function.
Maximize the area of polygon of N vertices and diameter of 1
The position of vertices is indicated as (rho,theta) coordinates
where rho denotes the distance to the base point
(vertex with number N) and theta the angle from the x-axis.
A user (black box) function is used to implement the problem.
(c) 2017 Fair Isaac Corporation
***********************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "xprs.h"
#include "xslp.h"
#define MAXROW 20
#define MAXCOL 20
#define MAXELT 50
#define MAXTOKEN 200
#define MAXCOEF 20
void XPRS_CC XSLPMessage(XSLPprob my_prob, void *my_object, char *msg,
int len, int msg_type);
double XPRS_CC MyFunc(double *Values, int *nArg);
int main(int argc, char *argv[]) {
XPRSprob xprob = NULL;
XSLPprob sprob = NULL;
int nRow, nCol, nSide, nElement, nToken, nCoef, nRowName, nColName;
int iRow, Sin, Cos;
char RowType[MAXROW];
double RHS[MAXROW], OBJ[MAXCOL], Element[MAXELT], Lower[MAXCOL], Upper[MAXCOL];
int ColStart[MAXCOL+1], RowIndex[MAXELT];
int ColIndex[MAXCOEF], FormulaStart[MAXCOEF+1];
int Type[MAXTOKEN];
double Value[MAXTOKEN], Factor[MAXCOEF];
int VarType[MAXCOL];
double InitialValue[MAXCOL];
int ReturnValue;
int i, j;
char RowNames[500], ColNames[500];
void *Func;
/* Initialisation */
if (ReturnValue=XPRSinit(NULL)) goto ErrorReturn;
if (ReturnValue=XSLPinit()) goto ErrorReturn;
if (ReturnValue=XPRScreateprob(&xprob)) goto ErrorReturn;
if (ReturnValue=XSLPcreateprob(&sprob,&xprob)) goto ErrorReturn;
/* XSLPsetcbmessage */
XSLPsetcbmessage(sprob,XSLPMessage,NULL);
nSide = 5;
nRowName = 0;
/* Rows */
nRow = nSide-2 + (nSide-1)*(nSide-2)/2 + 1;
for (i=0;i<nRow;i++) RHS[i] = 0;
nRow = 0;
RowType[nRow++] = 'E'; /* OBJEQ */
nRowName = nRowName + 1 + sprintf(&RowNames[nRowName],"OBJEQ");
for (i=1;i<nSide-1;i++) {
RowType[nRow++] = 'G'; /* T2T1 .. T4T3 */
RHS[i] = 0.001;
nRowName = nRowName + 1 + sprintf(&RowNames[nRowName],"T%dT%d",i+1,i);
}
for (i=1;i<nSide-1;i++) {
for (j=i+1;j<nSide;j++) {
RowType[nRow] = 'L';
RHS[nRow++] = 1.0;
nRowName = nRowName + 1 + sprintf(&RowNames[nRowName],"V%dV%d",i,j);
}
}
RowType[nRow] = '\0';
/* Columns */
nColName = 0;
nCol = (nSide-1)*2 + 2;
nElement = 0;
for (i=0;i<nCol;i++) {
OBJ[i] = 0; /* objective function */
Lower[i] = 0; /* lower bound normally zero */
Upper[i] = XPRS_PLUSINFINITY; /* upper bound infinity */
}
/* OBJX */
nCol = 0;
ColStart[nCol] = nElement;
OBJ[nCol] = 1.0;
Lower[nCol++] = XPRS_MINUSINFINITY; /* free column */
Element[nElement] = -1.0;
RowIndex[nElement++] = 0;
nColName = nColName + 1 + sprintf(&ColNames[nColName],"OBJX");
/* THETA1 - THETA 4 */
iRow = 0;
for (i=1;i<nSide;i++) {
nColName = nColName + 1 + sprintf(&ColNames[nColName],"THETA%d",i);
InitialValue[nCol] = 3.14159 * ((double) (i)) / ( (double) (nSide));
VarType[nCol] = 4;
ColStart[nCol++] = nElement;
if (i < nSide-1) {
Element[nElement] = -1;
RowIndex[nElement++] = iRow+1;
}
if (i > 1) {
Element[nElement] = 1;
RowIndex[nElement++] = iRow;
}
iRow++;
}
Upper[nCol-1] = 3.1415926;
/* Equals column */
nColName = nColName + 1 + sprintf(&ColNames[nColName],"=");
ColStart[nCol] = nElement;
Lower[nCol] = Upper[nCol] = 1.0; /* fixed at 1.0 */
InitialValue[nCol] = 1;
VarType[nCol] = 4 | 0x4000;
nCol++;
/* Remaining columns come later */
for (i=1;i<nSide;i++) {
Lower[nCol] = 0.01; /* lower bound */
Upper[nCol] = 1;
InitialValue[nCol] = 1;
VarType[nCol] = 4;
ColStart[nCol++] = nElement;
nColName = nColName + 1 + sprintf(&ColNames[nColName],"RHO%d",i);
}
ColStart[nCol] = nElement;
XPRSsetintcontrol(xprob,XPRS_CSTYLE,1);
XPRSsetintcontrol(xprob,XPRS_MPSNAMELENGTH,16);
if (ReturnValue=XPRSloadlp(xprob,"Polygon",nCol,nRow,RowType,RHS,NULL,OBJ,
ColStart,NULL,RowIndex,Element,Lower,Upper)) goto ErrorReturn;
if (ReturnValue=XPRSaddnames(xprob,1,RowNames,0,nRow-1)) goto ErrorReturn;
if (ReturnValue=XPRSaddnames(xprob,2,ColNames,0,nCol-1)) goto ErrorReturn;
/* Find index for SIN and COS */
if (ReturnValue=XSLPgetindex(sprob,XSLP_INTERNALFUNCNAMES,"SIN",&Sin)) goto ErrorReturn;
if (ReturnValue=XSLPgetindex(sprob,XSLP_INTERNALFUNCNAMES,"COS",&Cos)) goto ErrorReturn;
/* Define user function */
nToken = 0;
XSLPsetstring(sprob,&i,"MyFunc");
Type[nToken] = XSLP_STRING;
Value[nToken++] = (double) i;
Type[nToken] = XSLP_UFEXETYPE;
Value[nToken++] = (double) 0x01;
Type[nToken] = XSLP_UFARGTYPE;
Value[nToken++] = (double) 023;
XSLPsetstring(sprob,&i,"Dummy.DLL"); /* not important for internal function */
Type[nToken] = XSLP_STRING;
Value[nToken++] = (double) i;
Type[nToken] = XSLP_EOF;
XSLPloaduserfuncs(sprob,1,Type,Value);
XSLPaddnames(sprob,XSLP_USERFUNCNAMES,"MyArea",1,1);
Func = MyFunc;
XSLPchguserfuncaddress(sprob,1,&Func); /* redirect to internal function */
/* Build up nonlinear coefficients */
/* Area */
nToken = 0;
nCoef = 0;
RowIndex[nCoef] = 0;
ColIndex[nCoef] = nSide;
Factor[nCoef] = 1.0;
FormulaStart[nCoef++] = nToken;
Type[nToken] = XSLP_RB;
Value[nToken++] = 0;
for (i=nSide-1;i>0;i--) {
Type[nToken] = XSLP_COL;
Value[nToken++] = i;
Type[nToken] = XSLP_COL;
Value[nToken++] = nSide+i;
}
Type[nToken] = XSLP_FUN;
Value[nToken++] = 1;
Type[nToken] = XSLP_EOF;
Value[nToken++] = 0;
/* Distances */
for (i=1;i<nSide-1;i++) {
for (j=i+1;j<nSide;j++) {
RowIndex[nCoef] = iRow++;
ColIndex[nCoef] = nSide;
Factor[nCoef] = 1.0;
FormulaStart[nCoef++] = nToken;
Type[nToken] = XSLP_COL;
Value[nToken++] = nSide+i;
Type[nToken] = XSLP_CON;
Value[nToken++] = 2;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_EXPONENT;
Type[nToken] = XSLP_COL;
Value[nToken++] = nSide+j;
Type[nToken] = XSLP_CON;
Value[nToken++] = 2;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_EXPONENT;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_PLUS;
Type[nToken] = XSLP_CON;
Value[nToken++] = 2;
Type[nToken] = XSLP_COL;
Value[nToken++] = nSide+i;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_MULTIPLY;
Type[nToken] = XSLP_COL;
Value[nToken++] = nSide+j;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_MULTIPLY;
Type[nToken] = XSLP_RB;
Value[nToken++] = 0;
Type[nToken] = XSLP_COL;
Value[nToken++] = j;
Type[nToken] = XSLP_COL;
Value[nToken++] = i;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_MINUS;
Type[nToken] = XSLP_IFUN;
Value[nToken++] = Cos;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_MULTIPLY;
Type[nToken] = XSLP_OP;
Value[nToken++] = XSLP_MINUS;
Type[nToken] = XSLP_EOF;
Value[nToken++] = 0;
}
}
if (ReturnValue=XSLPloadcoefs(sprob,nCoef,RowIndex,ColIndex,Factor,
FormulaStart,1,Type,Value)) goto ErrorReturn;
for (i=0;i<nCol;i++) {
ColIndex[i] = i;
}
if (ReturnValue=XSLPloadvars(sprob,nCol-1,&ColIndex[1],&VarType[1],NULL,
NULL,NULL,&InitialValue[1],NULL)) goto ErrorReturn;
XSLPwriteprob(sprob,"Polygon_userfunc","");
if (ReturnValue=XSLPmaxim(sprob,"")) goto ErrorReturn;
if (ReturnValue=XSLPwriteslxsol(sprob, "Polygon_userfunc", "")) goto ErrorReturn;
goto NormalReturn;
ErrorReturn:
printf("\nError %d",ReturnValue);
NormalReturn:
XSLPdestroyprob(sprob);
XPRSdestroyprob(xprob);
XSLPfree();
XPRSfree();
return(ReturnValue);
}
void XPRS_CC XSLPMessage(XSLPprob my_prob, void *my_object, char *msg,
int len, int msg_type) {
switch (msg_type) {
case 4: /* error */
case 3: /* warning */
case 2: /* dialogue */
case 1: /* information */
printf("%s\n",msg);
break;
default: /* exiting */
fflush(stdout);
break;
}
}
double XPRS_CC MyFunc(double *Values, int *nArg) {
int i;
double Area;
Area = 0;
for (i=3;i<nArg[0];i=i+2) {
Area = Area + 0.5*Values[i-3]*Values[i-1]*sin(Values[i]-Values[i-2]);
}
return Area;
}
|
